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Summary 


We  begin  this  report  with  a  survey  of  direction-finding  algorithms  restricted  to  the  narrowband 
case  in  element  space.  We  then  set  out  to  search  for  a  simple  satisfactory  stress  measure  to 
direction-finding  algorithms  applied  to  scenarios  with  two  sources.  The  word  "stress”  for  a  scenario 
is  used  to  indicate  the  difficulty  that  an  algorithm  has  in  determining  the  direction  of  arrival 
of  each  of  the  two  sources.  A  scenario  is  said  to  be  "stressfu1”  to  an  algorithm  if  it  causes 
difficulty  in  resolving  the  sources  or  in  determining  their  directions  of  arrival.  The  apparent  lack 
of  understanding  of  the  effect  of  phase  difference  between  two  signals  in  direction-finding  made 
such  an  ideal  stress  measure  illusive.  Regardless  of  all  its  deficiencies,  eigenvalue  ratio  is  still  a 
useful  measure  in  summarizing  stressing  factors  in  many  aspects  for  two-arrival  direction- finding 
scenario. 

Eigenvalue  weighting  appears  in  some  noise  subspace  methods,  in  parametric  signal  subspace 
fitting  methods,  and  in  nonparametric  subspace  adaptive  nulling  beamforming.  For  a  two-source 
array-processing  scenario,  normalized  large  and  small  eigenvalues  A]  and  A2  are  reduced  to  forms 
depending  only  on  a  real  triplet:  phase-dependent  variable  £,  phase-independent  variable  77,  and 
power  ratio  The  pairs  (£,  77)  are  confined  to  an  isosceles-like  region.  We  characterize 

•  this  isosceks-like  region  and  the  many-to-one  mapping  from  the  Cartesian  product  of  the 
temporal  and  spatial  correlation  unit-disks  onto  this  region, 

•  the  behavior  of  the  eigenvalues  and  their  ratios  as  functions  of  the  real  triplet  both  analyt¬ 
ically  and  graphically  with  respect  to  direction-finding. 

The  main  contribution  of  this  work  is  a  manageable  presentation  of  a  compact  map  showing 
Aj,  A2,  and  ^  as  functions  of  £,  77,  and  ^  over  all  possible  scenarios.  This  enables  one  to  see 
the  relative  positions  among  different  scenarios.  We  also  present  some  easy- to- remember  formulas 
that  enable  one  to  exercise  ”back-of-envelope”  assessment  of  scenarios. 

The  small  eigenvalue  is  shown  to  diminish  qualitatively  and  quantitatively  for  two-arrival 
scenarios  increasingly  stressed  with  high  temporal  and/or  spatial  correlations.  The  special  case  of 
equal-strength  signal  arrivals  =  1),  also  important  in  low-angle  radar  tracking,  shares  many 
rich  structures  of  general  The  equal-strength  case  also  has  several  additional  unique  features 
for  signal  eigenvalue  ratio  j^-,  which  is  important  in  direction-finding: 

•  an  extra  6  dB  over  the  vertical  axis  £  =  0  as  compared  to  the  general  case, 

•  a  6-dB  increase  for  equipower  arrivals  highly  correlated  both  temporally  and  spatially  from 
changing  the  angle  difference  between  the  two  unit-disk  vectors  from  90°(orthogonal  wave¬ 
forms)  to  0°(in-phase  waveforms), 

•  the  lower  left  corner  of  the  isosceles-like  triangle  is  a  point  of  discontinuity  for  the  eigenvalue 
ratio  41- 

The  developed  results  are  used  to  assess  some  scenarios  used  by  Cadzow  and  Ottersten. 
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Chapter  1 

Introduction 


High-resolution  algorithms  for  direction-finding  (DF)  using  general  array  configurations  are  still 
maturing.  For  the  parametric  approach  to  the  narrowband  (NB)  problem  in  element  space,  some 
consensus  appears  to  be  emerging.  Therefore,  the  scope  of  this  report  on  direction-finding  is 
restricted  to  the  narrowband  case  in  element  space. 

Because  of  the  different  scenarios  used  among  the  reported  analytic  or  simulation  performance 
results,  the  evaluation  of  the  various  high-resolution  NB  DF  algorithms  is  made  difficult  by  the  lack 
of  a  common  footing  in  relating  and  assessing  them.  The  possible  combination  of  all  parameters  to 
specify  scenarios  to  test  and  evaluate  candidate  algorithms  is  just  too  vast.  Even  if  one  conducts 
such  an  exhaustive  evaluation,  the  effective  presentation  of  results  poses  another  problem.  The 
lack  of  any  simple  yardstick  to  measure  some  sort  of  scenario  stress  makes  it  difficult  to  effectively 
grasp  the  different  reported  results  and  to  form  a  broad  brush  mental  picture  summarizing  the 
very  complicated  scene.  The  word  ’’stress*  for  a  scenario  is  used  to  indicate  the  difficulty  that  an 
algorithm  has  in  determining  the  direction  of  arrival  of  each  of  the  two  sources.  A  scenario  is  said 
to  be  "stressful”  to  an  algorithm  if  it  causes  difficulty  in  resolving  the  sources  or  in  determining 
their  directions  of  arrival.  From  a  DF  user’s  point  of  view,  it  is  highly  desirable  to  have  some 
guidelines  deveioped  to  help  draw  an  envelope  of  scenario  parameters  corresponding  to  a  potential 
stressful  operating  environment. 

We  limit  the  scope  of  our  attempt  to  search  of  yardsticks  measuring  stress  of  DF  scenarios 
to  only  two  sources.  This  is  because  of  the  ease  of  obtaining  a  measure  to  evaluate  scenarios  by 
dealing  with  only  solutions  of  quadratic  characteristic  equation  and  that  liioo  t  Oi  the  reported  sim¬ 
ulation  results  in  the  literature  involve  only  two-source  scenarios.  Even  then  some  non-negligible 
effort  is  required  to  display  such  results. 

So  far,  we  found  such  a  stress  measure  illusive.  The  effect  on  direction-finding  of  the  phase 
difference  between  two  signals’  source  amplitudes  is  apparently  not  a  well-understood  subject. 
When  this  phase  difference  does  not  come  into  play  in  the  eigenvalue  expressions,  such  as  when 
the  two  arrivals  are  temporally  uncorrelated  or  spatially  orthogonal,  eigenvalue  ratio  seems  to  be 
a  plausible  stress  measure.  A  large  value  of  the  eigenvalue  ratio  will  reflect  the  cumulative  effect 
of  large  magnitudes  of  the  normalized  spatial  and  temporal  correlation  coefficients  and  the  power 
ratio  of  the  two  arrivals.  Most  reporting  on  the  phase  difference  effect  is  for  the  low-angle  radar 
tracking  problem  with  two  closely  spaced,  coherent,  equipowered  plane-wave  arrivals  impinging  on 
a  uniform  line  array.  In  examining  candidates  for  stress  measures,  we  have  considered  eigenvalue 
ratio,  Cramer-Rao  lower-bound  (CRLB),  and  asymptotic  performance  results  of  some  subspace 
algorithms  for  such  scenarios.  Somehow  we  cannot  reconcile  all  the  differences  among  these 
candidates  together  with  our  physical  intuition. 
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In  conducting  the  illusive  search  for  a  simple  satisfactory  stress  measure  for  scenarios  to 
direction-finding  algorithms,  we  start  the  main  body  of  chapter  1  with  a  survey  of  narrowband 
DF  algorithms  for  general  sensor  arrays  in  element  space.  This  includes  a  taxonomy  table,  char¬ 
acterization  and  asymptotic  performance  of  the  main  algorithm  categories,  and  the  initialization 
procedure,  which  is  the  critical  item  of  most  high-resolution  algorithms.  Then,  the  idea  of  for¬ 
mulating  guidelines  desired  for  DF  scenario  design  is  developed  for  the  following:  too  many 
parameters  for  a  two-source  DF  scenario  for  a  given  sensor  array,  simple  map(s)  to  show  relative 
positions  of  scenarios  and  yardstick  for  measuring  stress  of  scenarios  to  direction-finding  algo¬ 
rithms,  and  dominant  eigenvalues  chosen  as  the  key  items.  Next,  we  discuss  eigenvalues  for  the 
noise-free  covariance  matrix  of  two  discrete  sources  by  considering  the  following  issues:  treatment 
of  this  subject  in  Hudson's  text. [23] ,  the  original  6-dimensional  parameter  space  for  a  given  sensor 
array  under  the  noise-free  assumption,  compact  map(s)  another  manv-to-one  reduction,  and 
3-d  map  with  the  power  ratio  between  the  two  arrivals  as  third  coordinate. 

Chapter  2  addresses  the  notion  of  stress  to  direction-finding  algorithms  and  scenarios.  First, 
the  following  three  main  sources  for  the  signal  eigenvalue  spread  or  causes  for  small  signal  eigen¬ 
values  are  reviewed:  high  spatial  and/or  temporal  correlation  and/or  high  contrast  of  signal 
strength  between  the  two  arrivals.  Next,  the  eigenvalues  in  weighting  expressions  are  discussed 
for  three  categories:  noise  subspace  methods,  paiametric  signal  subspace  fitting  methods,  ami 
nonparametric  subspace  adaptive  nulling  beamforming.  Then,  the  drawbacks  of  the  oversimplify¬ 
ing  stress-measuring  yardstick  using  condition  numbers  are  mentioned:  infinite  condition  number 
for  all  100%-correlated  two  arrivals,  equal  condition  numbers  corresponding  to  different  CRLBs, 
CRLI3  not  appealing  to  intuitive  expectation,  and  the  redeeming  value  for  signal  eigenvalue  ratio. 
At  the  end,  we  explain  the  reasons  why  special  treatment  is  given  for  equal-power  arrivals. 

During  the  course  of  this  search  of  a  good  stress  measure,  we  thoroughly  characterized  the 
eigenvalues  and  their  ratio  for  two-source  scenarios  and  came  up  with  an  effective  way  to  present 
results.  Chapter  3  is  about  the  analytical  characterization  of  eigenvalues  for  two  arrivals.  Chap¬ 
ter  4  is  about  the  fractional  beainwidth  separation  between  two  plane  waves  impinging  on  a 
uniform  line  array  (ULA).  Chapter  5  is  about  the  plots  of  signal  eigenvalues  and  their  ratio  in  dll 
scab*.  The  large  and  small  eigenvalues  as  well  as  their  ratio  are  functions  of  a  phase-dependent 
variable  and  a  phase-independent  variable,  both  derived  from  the  normalized  spatial  and  temporal 
correlation  coefficients,  and  the  power  ratio  between  the  two  arrivals.  The  main  contribution  of 
this  report  is  a  manageable  presentation  of  a  compact  map  showing  these  three  functions  over  all 
possible  scenarios.  This  enables  one  to  see  the  relative  positions  among  different  scenarios.  YVe  also 
present  some  easy-to-remember  formulas  that  enable  us  to  exercise  back-of-enveiupe  assessment 
of  scenarios  during  presentations  or  discussions. 

These  developed  results  are  used  to  assess  some  scenarios  used  by  Cadzow  and  Ottersten  in 
chapter  6.  Chapter  7  serves  to  conclude  the  main  findings  of  this  the  report. 

1.1  A  Brief  Survey  of  Narrowband  DF  Algorithms  for  General 
Sensor  Arrays  in  Element  Space 

For  the  reader  interested  in  the  I)F  algorithms  for  general  sensor  arrays,  a  classification  chart  of 
some  of  the  candidate  narrowband  algorithms  is  helpful.  To  this  end.  we  provide  table  1.1.  Some 
brief  discussions  of  this  taxonomy  table  are  provided  in  the  following  subsections.  They  are  about 
characterization,  performance,  and  implementations. 
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1.1.1  A  Broad-Brush  Taxonomy  Table 

Over-simplification  and  oversight  is  unavoidable,  but  this  table  serves  as  a  starting  point.  Hope¬ 
fully,  the  discussion  will  remedy  some  of  the  simplistic  presentation.  The  reader  is  cautioned 
not  to  read  too  much  into  the  performance  comparison  parts  of  this  table.  A  single  measure 
of  resolution  performance,  sometimes  based  on  asymptotic  performance  expressions  and  often 
evaluated  subjectively,  using  symbols  +  or  -,  is  employed  here  to  compare  entries  in  the  same 
columns  or  rows.  This  appears  strange  as  the  report  discusses  a  multitude  of  scenarios  for  even 
a  two-arrival  case.  The  measure  is  even  iess  meaningful  when  entries  not  in  the  same  column  or 
row  are  compared.  A  broad-brush  table  such  as  this,  together  with  only  a  few  pages  of  discussion, 
definitely  will  not  suffice  for  some  workers  in  the  field.  The  complexity  involved  in  just  comparing 
the  simulation  work  from  two  researchers  can  be  inferred  by  reading  chapter  6,  which  is  about 
some  of  Cadzow’s  and  Ottersten’s  scenarios.  One  alternative  is  a  similar  large-scale  study  of  such 
scope  as  the  1984  MIT/Lincoln  Lab  report  [3]. 

1.1.2  Characterization  of  Some  Algorithms 

Ottersten  and  Viberg  [31,35,36,38,59,60.61]  jointly  formulated  a  unified  subspace  fitting  frame¬ 
work  that  relates  several  existing  methods  including  Conventional  Beam  Forming  (CBF),  Multiple 
Signal  Classification  (MUSIC),  multidimensional  MUSIC,  and  Deterministic  Maximum  Likelihood 
(I)et-ML).  Directions  of  arrivals  (DOAs)  are  estimated  individually  by  CBF  and  MUSIC;  i.e..  by 
searching  over  a  grid  in  a  nonparametric  fashion  as  is  done  by  the  Minimum  Variance  Distortion¬ 
less  Response  (MVDR)  Adaptive  Beamforming  (ABF)  algorithms.  DOAs  are  estimated  jointly  by 
multidimensional  MUSIC  and  Det-ML  by  using  parametric  methods  without  a  grid.  Those  using 
a  grid  cannot  deal  with  coherent  (100/L  correlation)  arrivals  for  general  array  configurations. 

The  generalized  MUSIC  by  Ferrara  and  Parks  (18]  for  diverse  polarizations,  and  by  Haber  and 
Zoltowski  [21]  for  coherent  multipath,  is  a  hybrid  of  the  above  two  categories  but  is  still  of  the 
spectrum  type;  i.e.,  the  search  is  over  a  product  DOA  parameter  space.  The  dimension  of  the 
product  space  depends  on  the  maximum  number  of  multipaths  modeled  for  the  emitters  within 
the  coverage.  For  each  grid  point  in  the  product  domain,  we  consider  all  linear  combinations  of 
known  multipaths  corresponding  to  the  product  angle  coordinates.  The  linear  combination  that 
has  the  minimum  projection  onto  the  noise  subspace  is  found.  Tim  reciprocal  of  the  magnitude- 
squared  of  this  quantity  is  defined  as  the  spectrum.  The  spectrum  at  each  grid  point  involves 
the  computation  of  the  smallest  generalized  eigenvalue  of  a  small  matrix,  the  size  of  which  is 
the  number  of  multipaths  modeled.  The  display  and  extraction  of  parameters  for  the  generalized 
MUSK’  are  difficult  when  the  dimension  of  the  underlying  parameter  space  for  the  product  domain 
exceeds  two,  i.e.,  except  for  the  simplest  azimuthal-only  scenarios  having  at  most  two  multipaths 
from  each  emitter. 

Even  though  we  try  to  limit  consideration  to  algorithms  applicable  to  general  array  configura¬ 
tions  and  signal  models,  an  exception  is  made  for  Estimation  of  Signal  Parameters  via  Rotational 
Invariance  Techniques  (ESPRIT)  [45,65].  This  is  because  of  the  many  attractive  features  of  ES¬ 
PRIT  for  doublet  array  configurations  including  the  case  of  plane  wave  arrivals  impinging  on 
a  ULA.  The  main  features  are  the  simplicity  of  computation  and  the  ability  to  provide  signal 
parameter  estimates  in  the  absence  of  detailed  knowledge  of  the  sensor  array  characteristics.  Be¬ 
cause  of  the  misconception  that  ESPRIT  avoids  the  calibration  problem  entirely,  Zoltowski  and 
Stavrinides  [65]  comment  that  the  required  calibration  involved  in  the  doublet  arrays  is  not  a 
trivial  task.  The  inherent  redundancy  built  into  the  ideal  ESPRIT  array  structure,  which  causes 
the  degraded  resolution  for  a  ULA  using  Root  MUSIC,  ran  compensate  for  imperfect  calibration 
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so  that  corresponding  sensors  in  the  doublet  arrays  need  be  only  approximately  identical  [65]. 

There  are  some  convergence  in  thoughts  about  the  MVDR,  and  the  subspace  or  MUSIC 
algorithms,  especially  in  the  sonar  community.  This  is  discussed  later  in  subsection  2.2.3. 

1.1.3  Asymptotic  Performance 

The  asymptotic  properties  of  the  parameter  estimates  from  a  general  weighted  fitting  criterion 
were  derived  by  Ottersten  and  Viberg,  and  by  applying  an  appropriate  weighting  and  parame¬ 
terization,  they  obtained  distributions  of  different  estimators.  They  minimized  the  asymptotic 
variance  with  respect  to  the  subspace  weighting,  leading  to  a  novel  method  termed  Weighted 
Subspace  Fitting  (WSF).  They  also  examined  the  asymptotic  behavior  of  the  Stochastic  Maxi¬ 
mum  Likelihood  (Sto-ML)  and  the  distribution  is  shown  to  be  the  same  as  for  WSF.  Therefore, 
WSF  achieves  the  best  resolution  performance  attainable  by  Sto-ML,  which  requires  a  nonlinear 
(cP  +  d-\-  l)-dim  optimization  where  d  is  the  number  of  modeled  arrivals  for  azimuthal  problems, 
with  the  much  smaller  burden  of  only  a  d-dim  one  similar  to  that  of  Det-ML.  For  uncorrelated 
arrivals  such  as  in  matched-field  processing  (MFP),  the  d2  term  is  reduced  to  d.  For  a  plane-wave 
model  with  both  azimuthal  and  elevation  measurements,  the  original  linear  term  d  becomes  2d. 

The  MUSIC  algorithm  is  a  ’’one  arrival  at  a  time”  search  of  the  parameter  space,  and  this 
inherently  leads  to  a  large  bias  due  to  the  finite  sample,  which  dominates  the  variance  until  the 
amount  of  data  is  very  large.  The  multidimensional  techniques  (those  using  joint  estimation  of 
multiple  arrivals)  have  a  smaller  bias  and  the  asymptotic  expressions  can  even  be  used  for  smaller 
samples. 

The  asymptotic  performance  analysis  by  Ottersten  and  Viberg  patterned  after  Stoica  and 
Nehorai  [53,54,55]  who  had  many  results  similar  to  those  of  Ottersten  and  Viberg.  The  Signal 
Eigenvector  (SE)  approach  by  Cadzow  [6,7]  uses  only  dominant  eigenvectors.  Therefore,  according 
to  Ottersten  (private  communication  with  the  author  of  this  report),  its  finite  sample  performance 
has  nonzero  bias  and  the  mean  square  error  may  have  the  advantage  of  those  of  the  reduced 
rank  principle  by  Scharf  and  Tufts  [47,48];  i.e.,  bias  is  traded  for  reduced  variance.  A  paper 
has  been  submitted  by  H.  Wang  and  G.H.  Wakefield,  on  nonasymptotic  performance  analysis  of 
eigenstructure-based  direction-of-arrival  estimation  using  perturbation  theory. 

1.1.4  Initialization;  The  Critical  Item  in  Implementation 

The  jointly  estimating  types  of  algorithms,  with  the  exception  of  Stochastic  Maximum  Likelihood 
Estimation  (Sto-MLE),  involve  fitting  a  linear  combination  of  steering  vectors  of  unknown  DOA 
parameters  to  data.  They  all  solve  a  separable  nonlinear  least  squares  problem  [20,26,46]  or 
its  equivalent,  sometimes  called  reduced  or  concentrated  form.  They  differ  only  in  the  right- 
hand  s'des  of  the  fitting  equations.  For  Deterministic  Maximum  Likelihood  Estimation  (Det- 
MLE),  the  right-hand  sides  are  snapshot  data  vectors.  For  SE  and  WSF,  the  right-hand  sides 
are  the  unweighted  and  weighted  signal  eigenvectors  respectively.  They  can  differ  in  performance 
considerably.  Computationally,  these  approaches  are  similar  to  each  other  in  that  some  Gauss- 
Newton  algorithms,  or  their  approximations,  can  solve  the  parameter  fitting  problem  very  well  if  a 
good  initial  estimate  is  available;  i.e.,  one  close  to  one  of  the  global  minimum  (nonunique  to  within 
permutations).  Note  that  the  DOA  parameters  for  the  multiple  emitters  are  estimated  jointly  as 
they  correspond  to  a  single  point  in  the  parameter  space  over  which  the  least  squares  problem 
is  carried  out.  This  is  contrasted  to  the  spectrum  estimation  methods  where  each  emitter’s  (or 
arrival’s)  parameter  values  are  read  off  one  at  a  time. 


Getting  a  good  initial  estimate  is  definitely  no  small  feat,  as  the  cost  function  of  the  global 
minimum  search  has  many  relative  or  local  minima.  The  abstract  by  Degerine  as  announced  in 
the  manuscripts  in  review  of  the  August  1990  IEEE  ASSP  Transaction  [14]  has  the  following: 
The  likelihood  function,  in  estimating  a  Toeplitz  covariance  matrix,  can  present  a  local  maximum 
close  to  the  global  maximum.  An  analytical  proof  is  conducted  for  the  matrices  of  order  3  and 
numerical  examples  are  given  for  orders  3,  4,  and  5. 

Both  Ottersten  and  Viberg,  and  Cadzow  used  the  initialization  scheme  by  Ziskind  and  Wax  [63] 
to  tackle  this  critical  problem.  The  major  share  of  the  overall  direction-finding  computation  is 
essentially  from  this  initialization.  Starting  from  a  CBF  to  choose  the  first  strongest  arrival,  the 
Gram-Schmidt  orthogonal  procedure  is  successively  applied  to  form  (multiple)  null(s)  at  previously 
selected  arrival(s)  [23,  pages  41-46]  to  choose  the  next  strongest  arrival.  In  a  two-signal  scenario 
we  sized  for  azimuthal  and  elevation  angle  measurements  for  a  Mill-cross  array,  the  initial  estimate 
is  99%  of  the  total  computation  count. 

There  have  been  published  works  promoting  the  use  of  Expectation- Maximization  (EM)  algo¬ 
rithm  or  its  variants  for  the  direction-finding  problem.  While  this  was  met  with  skepticism  during 
off-the-record  discussions,  as  usual  such  skepticism  was  not  openly  reported.  But  the  recent  dis¬ 
sertation  by  Wu  [62]  shows,  through  both  simulation  and  theoretical  analysis,  that  conventional 
methods  such  as  the  gradient  and  the  Newton  methods  have  more  desirable  performance  in  (local) 
convergence  than  the  EM  algorithm. 

Note  that  the  EM  algorithm  and  its  variant  depend  equally  on  the  choice  of  a  good  initial 
estimate  as  do  the  more  conventional  methods.  In  fact,  some  of  these  EM  numerical  experiments 
such  as  in  Wu’s  dissertation  [62]  used  exactly  the  same  Ziskind  and  Wax’s  algorithm  to  find  an 
acceptable  initial  estimate.  In  a  recent  paper,  Ziskind  and  Wax  [64]  incorporated  some  a  priori 
periodicity  structural  information  of  the  spatial  covariance  function  for  a  ULA  with  uncorrelated 
arrivals  in  a  Det-ML  scheme.  They  obtained  higher  resolution  than  MUSIC  but  used  these 
’’lower”  resolution  algorithms  to  initialize  the  EM  method.  (Interestingly,  Ziskind  and  Wax  did 
not  use  the  initialization  scheme  bearing  their  own  names.)  Lo,  Nagaraj,  and  Rukhin  [31]  also 
used  the  Ziskind  and  Wax  [63]  initialization  scheme  in  their  DF  algorithm  using  cyclic  regression 
for  WSF.  While  their  claim  of  great  reduction  of  computational  complexity  is  highly  debatable, 
they  obtained  interesting  results  by  modifying  Det-ML  to  achieve  an  asymptotically  equivalent 
expression  of  the  WSF  criterion.  We  are  somewhat  sold  on  the  principle  of  rank-reduction  and 
still  consider  the  rank-reduced  version  of  the  original  WSF  as  the  favored  baseline  candidate. 

As  the  initialization  is  where  resource  should  be  allocated  and  the  description  of  the  most 
commonly  used  Ziskind  and  Wax’s  algorithm  was  terse,  we  characterize  it  here  to  provide  insights. 
Hopefully,  someone  can  provide  less  expensive  alternatives  in  the  future.  Starting  from  a  CBF  to 
choose  the  first  strongest  arrival,  the  initialization  continues  in  a  recursive  fashion  using  successive 
adaptive  nulling.  Suppose  k-1  strongest  arrivals  have  been  selected.  The  next  or  the  k-th  strongest 
arrival  is  selected  from  the  highest  peak  of  the  next  adaptive  nulling  spectrum  in  the  following  way. 
The  weight  vector  of  the  adaptive  nulling  at  each  grid  point  is  made  orthogonal  to  the  steering 
vectors  of  all  the  k-1  strong  arrivals  previously  selected  to  produce  a  distortionless  response  for 
an  arrival  with  the  DOA  parameters  corresponding  to  this  grid  point.  Because  of  the  adaptive 
nulling  application,  the  resolution  of  this  initialization  is  at  least  medium  high,  so  a  grid  finer  than 
that  for  the  CBF  is  needed.  Such  orthogonal  beamforming  effectively  subtracts  out  strong  signals 
from  sensor  data  by  using  orthogonal  projection  without  performing  the  dangerous  subtraction 
of  large  numbers.  Ziskind  and  Wax  used  a  projection  matrix  update  formula  based  on  sequential 
one-step  application  of  the  Gram-Schmidt  procedure.  As  the  simplicity  of  this  algorithm  derives 
from  initializing  or  updating  only  one  source’s  DOA  parameters  at  each  adaptive  nulling  scan,  it 
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is  similar  to  the  coordinate  descent  method  well-known  in  multivariate  optimization.  The  original 
paper  by  Ziskind  and  Wax  went  beyond  just  initialization.  The  same  procedure  was  proposed  to 
continue  cyclically,  i.e.,  wrapping  around  the  source  indices:  1, if  d  contacts  are  modeled. 
This  DF  algorithm  proper  has  not  been  received  well  because  of  the  exceedingly  slow  convergence 
rate,  which  is  no  surprise  given  the  reputation  of  coordinate  descent  methods.  However,  the 
similarity  between  the  coordinate  descent  and  Ziskind  and  Wax’s  methods  stops  here,  as  Ziskind 
and  Wax’s  initialization  aims  at  a  global  multivariate  search  while  coordinate  descent  method 
does  not. 

Some  potential  improvements  on  the  Ziskind  and  Wax  initialization  scheme  are  obvious.  For 
those  sources  interfering  with  each  other  at  most  marginally  because  of  wide  spatial  separations 
and  the  lack  of  strength  contrast,  one  may  entertain  the  idea  of  selecting  them  simultaneously  in 
one  adaptive  nulling  scan.  On  the  other  hand,  for  a  pair  of  nearby  sources  interfering  with  each 
other  within  a  mainlobe,  a  fixed  grid  spacing  throughout  the  sequential  orthogonal  beamforming 
may  not  serve  the  purpose  well.  We  may  refine  the  grid  spacing  near  the  selected  strong  sources 
so  that  potentially  present  nearby  sources  may  be  selected.  We  may  also  perturb  the  grid  system 
locally  by  offsetting  some  grid  lines  because  we  cannot  reuse  for  later  adaptive  nulling  scans  those 
old  grid  points  already  occupied  by  selected  strong  sources. 

The  delineation  shown  in  the  taxonomy  chart  between  subspace  DF  and  ABF  is  of  broad-brush 
nature.  For  MUSIC-type  algorithms,  there  are  proposed  implementations  without  the  need  to 
perform  eigen-decompositions.  The  works  by  Munier  and  Delisle  [33],  Reilly,  Chen,  and  Wong  [43], 
and  Friedlander  [19]  all  appeared  in  1988  within  a  short  duration.  On  the  other  hand,  the  traffic 
is  not  just  one-way.  Clark  and  Roberts  [12]  advocated  using  rank-one  eigenstructure  updating  in 
realtime  ABF. 


1.2  A  Yardstick  to  Measure  Stress  of  a  Two-source  Scenario 
Exacting  on  a  DF  Algorithm 

1.2.1  Guideline  Desired  for  DF  Scenario  Design 

From  the  perspectives  of  High-Frequency  Direction-Finding  (HFDF)  system  development  and  the 
user’s  point  of  view,  it  is  important  to  test  and  evaluate  the  various  algorithms  against  realistic 
scenarios  by  using  either  measured  or  simulated  data.  These  scenarios  can  be  used  to  evaluate 
the  relevancy  and  deficiency  of  the  basic  assumptions  underlying  algorithm  development,  and 
to  calibrate  simulation  results.  The  cost  of  collecting  measured  data  determines  the  number  of 
scenarios  possible.  Therefore,  guidelines  are  needed  to  select  scenarios  to  meet  the  time  and 
budget  constraints  for  data  collection  and  postanalysis. 

1.2.2  Too  Many  Parameters  for  Two-Source  DF  Scenarios  for  a  Given  Sensor 
Array 

Even  for  highly  idealized  scenarios  with  just  two  discrete  sources  and  a  given  array  with  perfectly 
specified  calibration,  there  are  too  many  parameters  to  be  exhaustively  covered.  These  parameters 
are  the  two  signal-to-noise  ratios  (SNRs),  number  of  snapshots,  source  correlation  coefficient,  and 
DOAs  of  sources. 

As  a  real  data  set  is  difficult  to  obtain  and  the  understanding  of  algorithms  axe  still  maturing, 
a  performance  study  using  simulation  is  the  only  vehicle  available  for  most  algorithm  investiga¬ 
tors.  For  each  chosen  set  of  the  above  mentioned  parameters,  a  large  Monte  Carlo  simulation  is 
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needed  to  evaluate  each  candidate  algorithms  with  their  various  algorithm  settings.  Analytical 
performance  predictions,  when  available,  can  significantly  reduce  the  cost  of  such  Monte  Carlo 
simulations  and  sometimes  provide  understanding  of  the  algorithms  not  possible  to  obtain  by 
simulation.  Even  then,  the  number  of  possible  scenarios  makes  comparison  of  algorithms  difficult. 

1.2.3  Simple  Map(s)  to  Show  Relative  Positions  of  Scenarios  and  Yardsticks 
for  Measuring  Stress  of  Scenarios  to  Direction-Finding  Algorithms 

Typically,  each  author  of  a  particular  algorithm  tends  to  present  scenarios  favorable  to  the  al¬ 
gorithms  he  or  she  (re)invented.  Objectivity  and  modesty  is  not  pervasive.  Even  for  the  many 
well-conducted  comparisons  separately  conducted,  it  is  difficult  to  relate  these  seemingly  "apple 
and  orange”  scenarios.  Because  of  the  bodies  of  reported  results  available  in  the  literature  and 
the  cost  of  duplicating  such  tests  for  even  limited  standardized  scenarios,  we  settle  for  a  com¬ 
bination  of  simple  map(s)  and  a  realistic  yardstick.  We  mentioned  earlier  that  there  are  too 
many  parameters  to  be  covered  for  even  a  two-source  scenario.  The  compact  map  or  maps  will 
enable  us  to  show,  on  some  common  footing,  the  relative  positions  of  the  noise-free  version  of 
the  various  scenarios  already  reported  or  yet  to  be  used.  The  yardstick  is  intended  for  measuring 
the  stress  factor  of  the  scenarios  presented  to  the  DF  algorithms  already  reported  or  two-source 
simple  scenarios  to  be  planned  in  the  future.  Such  a  measure  does  not  account  for  the  cost  issues, 
such  as  operation  count  or  implementation  ease,  which  are  very  important  algorithmic  issues  in  a 
DF  system.  But  it  connects  to  the  issues  related  to  resolution  or  standard  deviation  of  the  DOA 
estimates  when  the  noise  floor  is  considered.  The  relative  level  of  the  noise  floor  is  determined  by 
the  SNR  and  the  number  of  snapshots  available. 

Hopefully,  we  can  relate  a  sonar  scenario  using  a  3-d  random  array  with  azimuth  and  eleva¬ 
tion  angle  measurements  reported  in  one  publication  to  a  HF  scenario  using  U LA  with  azimuth 
measurement  only  reported  in  another  publication.  We  leave  the  actual  performance  evaluation 
of  the  various  candidate  algorithms  to  either  analytical  or  numerical  methods  reported  or  yet  to 
be  carried  out.  What  we  provide  is  a  common  footing  for  comparing  the  noise-free  versions  of 
one  scenario  to  another  using  a  stress  measure. 

1.2.4  Dominant  Eigenvalues  Chosen  as  the  Key  Items 

The  key  components  of  a  DF  problem  are  the  detection  and  estimation  of  dominant  arrivals. 
Evidence  appeared  recently  that  the  WSF  approach  is  the  leading  candidate  for  a  baseline  high- 
resolution  HFDF  algorithm  for  a  general  sensor  array.  Furthermore,  Viberg  and  Ottersten  have 
provided  a  unified  frame  for  the  many  subspace  and  Det-MLE  algorithms. 

The  various  algorithms  are  different  in  this  unified  frame  in  the  weights  used.  These  weights 
are  expressions  of  the  eigenvalues  of  the  data  covariance  matrices.  To  this  end,  we  investigated 
the  eigenvalues  for  noise-free  two-source  problems. 


1.3  Eigenvalues  for  the  Noise- Free  Covariance  Matrix  of  Two 
Discrete  Sources 

1.3.1  Treatment  of  this  Subject  in  Hudson’s  Text 

Hudson’s  text  on  ABF  contains  a  four-page  section[23,  pages  52-55]  on  the  eigenvalues  for  the  data 
covariance  matrix  of  two  discrete  sources.  He  provides  the  general  formula  of  the  two  eigenvalues 


8 


and  some  special  case  treatment.  Apparently,  eigen-decomposition  methods  were  not  mature 
then,  the  treatment  did  not  realize  its  full  potential.  Here  we  interpret  Hudson’s  expressions 
differently  and  use  graphs  from  our  yardstick’s  point  of  view  to  provide  compact  maps  to  relate 
the  various  scenarios  and  to  provide  a  ”stress”  measure  that  can  ultimately  cause  some  algorithms 
to  break  down  in  the  standard  deviation  performance  of  their  DOA  estimates.  We  adhere  as  much 
as  possible  to  Hudson’s  original  notations  in  his  four-page  treatment. 

1.3.2  The  Original  Six-Dimensional  Parameter  Space  for  a  Given  Sensor  Array 
under  the  Noise-Free  Assumption 

For  a  given  sensor  array  with  noise-free  assumption  for  two  sources,  the  eigenvalues  depend  upon 
the  following  six  parameters: 

•  xi  and  7T2,  signal  powers  for  signals  1  and  2  with  tt i  >  x 2 

•  />,  normalized  multipath  correlation  coefficient  between  the  two  signals’  complex  amplitude, 
a  complex  number  constrained  on  a  unit-disk, 

•  <j>,  normalized  spatial  correlation  coefficient  between  the  two  arrivals’  steering  vectors,  also 
a  complex  number  constrained  on  a  unit-disk. 

For  simplicity,  we  may  assume  that  the  eigenvalues  have  been  scaled  with  respect  to  the 
product  of  the  number  of  sensors  and  7^.  This  way,  we  reduce  the  number  of  parameters  to 
five;  i.e.,  ^  and  an  element  from  a  four-dimensional  region  formed  by  the  Cartesian  product  of 
two  unit-disks.  We  have  achieved  a  many-to-one  mapping  from  the  set  of  all  possible  noise-free 
scenarios  to  the  admissible  set  of  five  parameters. 

While  we  let  the  point  on  each  unit-disk  be  arbitrary,  this  may  not  be  the  case  for  specific 
scenarios.  In  the  case  of  the  spatial  correlation  disk  for  a  ULA  with  a  given  number  of  elements, 
the  normalized  complex  inner  product  of  two  steering  vectors  is  a  (real)  sine  function  with  an 
associated  phase  term.  The  modulus  and  argument  are  related  in  some  way  and  together  trace 
out  a  curve  on  the  polar  plot  on  the  unit-disk.  Such  curves  depend  upon  the  number  of  sensors 
used.  Similar  cases  may  be  made  on  the  multipath  correlation  unit-disk,  perhaps  by  using  some 
argument  from  ray-tracing  or  other  propagation  models.  Also,  we  will  see  that  the  stress  factor 
is  not  a  strong  function  of  the  phase  terms  when  the  strength  ratio  is  away  from  unity. 

While  this  mapping  from  an  arbitrary  sensor  scenario  to  a  five- dimensional  data  point  in  the 
parameter  space  accomplishes  a  significant  reduction,  it  is  still  not  easy  to  gain  insight  of  the 
problem.  Hudson  considered  some  important  special  cases,  but  several  important  features  were 
left  out. 

1.3.3  Compact  Map(s):  Another  Many-to-One  Reduction 

In  the  following,  we  provide  an  additional  many-to-one  mapping  such  that  for  a  given  value  of  the 
strength  ratio  the  Cartesian  product  of  the  two  unit-disks  can  be  mapped  into  an  isosceles-like 
region  in  two  new  parameters.  The  two  new  parameters  are: 

(.  =  \p\\4>\  cos(argp  -  arg<£), 

=  M2)(i  -  M2)- 

For  convenience,  we  call  t?  the  phase-independent  variable  and  £  the  phase- depen  dent  variable. 
The  rationale  for  the  names  is  provided  in  section  3.1. 
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Note  the  symmetry  between  p  and  <f>  in  the  above  expressions.  Given  ir\  and  7T2,  they  are  the 
only  expressions  needed  in  computing  the  general  formula  for  the  two  eigenvalues.  These  two  new 
parameters  are  not  entirely  independent  of  each  other  as  will  be  explained  in  detail  later.  This  is 
in  contrast  to  the  original  six  parameters,  which  could  be  varied  independently  with  some  minor 
qualifying  remarks. 

We  can  further  identify  points  corresponding  to  the  following  four  special  cases  that  were 
pointed  out  by  Hudson: 

•  Arrivals  0%  correlated  with  respect  to  both  multipath  and  steering  vectors, 

•  Arrivals  0%  correlated  with  respect  to  either  multipath  or  steering  vectors  but  not  to  both, 

•  Arrivals  0%  correlated  with  respect  to  either  multipath  or  steering  vectors  and  100%  corre¬ 
lated  with  respect  to  the  other, 

•  Constructive  interference  with  the  two  arrivals  100%  correlated  with  respect  to  both  multi- 
path  and  steering  vectors  with  cosine  term  +1. 

Because  of  the  nature  of  this  many-to-one  mapping,  the  same  (£,rj)  point  can  come  from  several 
different  special  cases  in  the  original  four-dimensional  parameter  region  that  is  a  Cartesian  product 
of  two  unit-disks.  For  example,  the  point  (£,t?)  =  (0,0)  corresponding  to  the  special  case: 

Arrivals  0%  correlated  with  respect  to  either  multipath  or  steering  vectors  and  100% 
correlated  with  respect  to  the  other, 

* 

also  corresponds  to  the  following  special  case: 

’’Quadrature  summed”  interference  with  the  two  arrivals  100%  correlated  with  respect 
to  both  multipath  and  steering  vectors  with  zero  phase-difference  term;  i.e.,  the  two 
multipath  sources  are  in  the  same  direction  but  are  90°  out  of  phase. 

The  vertical  symmetry  line,  £  =  0,  corresponding  to  the  special  cases: 

Arrivals  0%  correlated  with  respect  to  either  multipath  or  steering  vectors  but  not  to 
both 

also  corresponds  to  the  following  special  cases: 

Arrivals  with  the  phasor  directions  on  the  two  unit-disks  perpendicular  to  each  other; 
i.e.,  argp  -  argd>  =  90°. 

Hudson  also  ignored  the  important  special  case: 

Destructive  interference  with  the  two  arrivals  100%  correlated  with  respect  to  both 
multipath  and  steering  vectors  and  with  the  cosine  term  -1. 

The  whole  baseline  of  the  isosceles-like  region  can  be  interpreted  as  special  cases  with  the 
two  arrivals  100%  correlated  with  respect  to  either  multipath  or  steering  vectors.  It  can  also 
be  interpreted  as  special  cases  with  the  two  arrivals  100%  correlated  with  respect  to  both  the 
multipaths  and  steering  vectors  in  addition  to  having  a  proper  valued  cosine  term.  We  note  that 
Hudson  unnecessarily  restricted  this  just  mentioned  special  case  to  the  single  point  corresponding 
to  the  constructive  interference. 
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We  also  provide  the  characterizations  with  respect  to  the  cosine  term  in  the  ^-definition.  The 
right-half  of  the  isosceles  region  corresponds  to  arrivals  with  the  phase  on  the  two  unit-disks 
differing  less  than  90°.  The  left-half  of  the  isosceles  region  corresponds  to  arrivals  with  the  phase 
on  the  two  unit-disks  differing  more  than  90°. 

Later  in  figure  3.2,  we  will  see  that  Hudson’s  textbook  presented  4  out  of  the  17  characteriza¬ 
tions  shown  there.  These  4  cases  are  denoted  by  asterisks. 

1.3.4  A  3-D  Map  with  Power  Ratio  ^  as  the  Third  Coordinate? 

We  can  add  the  ratio  as  the  third  coordinate.  However,  the  observations  made  after  examining 
the  three  surfaces  of  Aj,  A2,  and  their  ratio  indicate  that  their  dependence  of  ^  on  the  ratio  ^  is 
essentially  translational  in  the  dB  scale  especially  away  from  the  unity  ratio  base  plane.  Therefore, 
the  treatment  of  the  parameter  ^  should  be  handled  more  economically  as  an  additive  term  ^  in 
the  dB  scale.  In  other  words,  there  is  no  need  to  go  to  a  three-dimensional  map,  if  ever  possible. 
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Chapter  2 

In  Search  of  a  Simple  Stress 
Measure  for  Direction-Finding 
Algorithms  and  Scenarios 


We  will  first  discuss  the  three  main  sources  for  the  signal  eigenvalue  spread  or  causes  for  small 
signal  eigenvalues.  Next,  we  elaborate  on  eigenvector  weighting  explosions:  in  noise  subspace 
methods,  in  parametric  signal  subspace  fitting  methods,  and  in  nonparametric  subspace  adap¬ 
tive  beamforming.  Finally  we  mention  the  drawbacks  of  oversimplifying  the  stress  measuring 
yardstick:  infinite  condition  number  for  all  100%-correlated  two  arrivals,  equal  condition  num¬ 
bers  corresponding  to  different  Cramer-Rao  lower  bounds  (CRLBs),  and  CRLB  not  appealing  to 
intuitive  expectation. 


2.1  The  Three  Main  Sources  for  the  Signal  Eigenvalue  Spread 
or  Causes  for  Small  Signal  Eigenvalues 

Schmidt’s  dissertation  qualitatively  pointed  out  the  following  three  stressing  factors  to  a  direction¬ 
finding  algorithm.  For  small-  or  medium-sized  sensor  arrays,  the  resolution  is  not  sufficiently  high 
so  we  can  have  unresolved  arrivals.  Furthermore  100%  source  correlation  will  cause  rank  deficiency 
of  the  noise-free  data  covariance  matrix  as  well.  The  other  stress  factor  to  a  direction-finding 
algorithm  or  scenario  is  the  strength  ratio  between  the  two  signals. 

Hudson  presented  some  eigenvalue  results[23,  pages  52-55]  of  the  sensor  data  covariance  matrix 
for  two  discrete  sources.  We  can  use  his  results  here  to  infer  that  it  is  not  unusual  to  find  a  30- 
dB  spread  in  the  signal  eigenvalues.  Much  more  detailed  treatment  will  be  provided  later  in 
this  report.  The  three  multiplicative  factors  that  contribute  to  the  eigenvalue  spread  are:  the 
signal-strength  ratio,  the  angular  proximity,  and  the  high  correlation  between  the  two  sources. 
A  high  signal-strength  ratio,  say  10  dB,  almost  directly  translates  into  the  eigenvalue  ratio.  The 
generalized  cosine  between  the  two  steering  vectors  corresponding  to  the  two  arrivals  is  equivalent 
to  the  source  correlation  coefficient,  as  far  as  their  effect  on  the  eigenvalue  spread  is  concerned.  A 
value  of  95%  source  correlation  will  contribute  to  10  dB  in  the  spread.  A  ULA  with  10  elements  at 
half-wavelength  spacing  and  a  2-degree  source  spacing  between  a  second  source  and  a  broadside 
one  will  yield  0.95  value  of  the  generalized  cosine  between  the  two  steering  vectors,  and  therefore, 
another  10-dB  spread. 

The  condition  number,  which  in  this  case  is  the  ratio  of  the  two  dominant  eigenvalues  of 
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the  sensor  data  covariance  matrix,  is  a  good  candidate  to  quantitatively  describe  the  cumulative 
stress  of  scenarios.  It  is  of  interest  to  know  how  the  eigenvalues  behave  when  spatially  resolved 
or  unresolved  (with  respect  to  conventional  beamforming)  and/or  temporally  highly  (including 
100%)  correlated  arrivals  are  present,  and  as  the  ratio  between  the  two  signals’  strength  increases. 


2.2  Eigenvector  Weighting  Expressions 

2.2.1  In  Noise  Subspace  Methods 

Johnson  and  DeGraaf  [25]  proposed  to  modify  the  original  noise  subspace  MUSIC  algorithm  by 
weighting  the  noise  eigenvectors  with  the  noise  eigenvalues.  The  modification  was  derived  in  a 
somewhat  ad  hoc  fashion.  Some  real  data  experience  trading  off  Johnson  and  DeGraaf’s  weighting 
scheme  to  Schmidt’s  unweighted  MUSIC  was  reported  by  Martin  [32].  Martin’s  result,  based  on 
measured  data,  contradicted  Johnson  and  DeGraaf’s  result,  which  was  based  on  simulation  and 
showed  promise  using  noise  eigenvector  weighting. 

As  Martin  did  not  give  an  explanation  of  why  this  was  so  from  using  a  specific  HFDF  data  set, 
we  make  the  following  conjecture.  The  noise  eigenvalues  appeared  in  the  cost  function  through 
a  setting  of  ’’reciprocal  of  reciprocals”  of  noise  eigenvalues.  When  dealing  with  real  data,  as 
contrasted  with  simulated  data  under  the  ideal  noise  eigenvalues  assumption,  those  small  and 
furthermore  nonequal  eigenvalues  may  indeed  play  havoc. 

On  the  other  hand,  this  specific  HFDF  data  set  was  severely  marred  by  nonrandom  calibration 
error.  Furthermore,  the  input  data  to  these  two  noise-subspace  processors  were  the  output  from 
a  de-chirping  process  because  the  signal  waveform  was  Linear  Frequency  Modulated  (LFM). 

What  could  be  cast  as  uncertain  or  inconclusive  result  from  studying  this  specific  HFDF  data 
set  is,  however,  remedied  by  the  recent  analytical  result  from  Stoica  and  Nehorai  [54]  which  shows 
that  the  minimum  variance  in  the  class  of  weighted  noise  subspace  MUSIC  estimators  is  achieved 
by  the  nonweighted  MUSIC. 

2.2.2  In  Parametric  Signal  Subspace  Fitting  Methods 

After  Schmidt’s  MUSIC  algorithm  [49],  Viberg  and  Ottersten’s  recent  work  [34,35,36,38,59,60,61] 
on  weighted  subspace  fitting  (WSF)  direction-finding  is  among  the  most  important  works  in 
array  processing.  In  earlier  multidimensional  MUSIC  algorithms,  the  residuals  from  fitting  each 
signal  eigenvectors  by  a  linear  combination  of  steering  vectors  with  unknown  DOA  parameters 
are  weighted  equally  in  the  separable  nonlinear  least  squares  formulation  to  jointly  estimate  these 
DOA  parameters.  Viberg  and  Ottersten  found  that  a  suitable  weight  for  each  eigenvector’s  fitting 
can  be  used  to  minimize  the  asymptotic  estimation  error  of  DOA  parameters. 

The  weighting  mentioned  here  is  different  from  Johnson  and  DeGraaf’s  noise  eigenvector 
weighting  [25],  because  here  signal  eigenvector  weighting  is  involved.  The  optimal  WSF  method  for 
parametric  array  processing  has  (dominant)  signal  eigenvalues  appearing  in  the  weight  expressions 
in  the  cost  function.  WSF  is  an  optimal  solution  using  a  well-defined  criterion.  We  note  that  signal 
eigenvalues  are  much  more  stable  than  the  noise  eigenvalues  both  numerically  and  statistically. 

For  high  SNRs,  the  weights  used  in  the  Det-ML  and  the  optimal  WSF  methods  are  very 
nearly  equal.  However,  Cadzow’s  simulation  results  [8]  using  two  coherent  sources  with  high 
SNRs  indicate  that  the  Det-ML,  which  should  be  close  to  WSF  in  performance,  is  outperformed 
by  the  single  signal  eigenvector  fitting  method.  During  a  private  communication  with  us  discussing 
this  apparent  discrepancy,  Ottersten  attributed  this  to  the  rank  reduction  principle  by  Scharf  and 
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Tufts  [47]  and  Scharf  [48].  The  rank  reduction  principle  trades  off  model  bias  and  model  variance 
in  the  analysis  and  synthesis  of  signals  to  minimize  the  mean  squared  error.  Even  though  the 
rank  reduced  approach  provides  biased  estimates,  the  mean  square  error  could  be  lower  than  the 
CRLB,  which  is  asymptotically  approached  by  Det-ML,  an  asymptotically  unbiased  estimator. 
This  suggests  that  only  dominant  eigenvectors  should  be  used. 

The  paper  by  Swindlehurst  and  Kailath  [57]  discusses  phase  and  amplitude  gain  uncertainty 
and  sensor  position  errors.  The  optimal  weights  are  different  from  the  ideal  calibration  case. 

2.2.3  In  Nonparametric  Subspace  Adaptive  Nulling  Beamforming 

We  need  a  brief  description  of  the  evolution  of  the  well  known  Sampled  Matrix  Inverse  (SMI) 
method  in  adaptive  beamforming  (ABF)  to  the  subspace  adaptive  nulling  method  before  we  can 
address  the  equivalent  version  of  weighted  subspace  variants. 

For  sonar  applications,  large  arrays  lead  to  problems  in  collecting  enough  snapshots  for  the 
SMI  method  within  the  constraint  of  stationarity.  Snapshot  sample  requirements  results  are 
expressed  in  terms  of  the  SNR  degradation  because  of  the  inaccurate  covariance  matrix  estimation. 
For  SMI,  the  number  of  snapshots  must  exceed  twice  the  number  of  sensors  if  the  tolerable 
SNR  degradation  is  to  be  no  more  than  3  dB.  Perhaps  SNR  is  better  written  as  SBR  (Sienal- 
to- Background  Ratio)  or  SINR  (Signal-to-Interference-and-Noise  ratio)  as  some  people  prefer, 
win  this  case,  the  increased  background  level  is  due  to  the  interference  left  over  because  of  the 
inadequate  adaptive  cancellation,  whether  in  the  presence  or  absence  of  signal.  The  net  result  is 
a  decrease  in  SNR  or  SBR.  which  in  turn  leads  to  the  system  performance  loss  such  as  probability 
of  detection  from  using  a  higher  absolute  detection  threshold  (with  respect  to  the  white  noise, 
e.g.)  to  maintain  a  fixed  probability  of  false  alarm. 

Subspace  adaptive  nulling  depends  upon  the  number  of  dominant  eigenvectors.  In  their  Prin¬ 
cipal  Component  Inverse  (PCI)  papers,  Kirsteins  and  Tufts  [28,29]  concluded  that  the  number  of 
snapshots  needs  to  exceed  only  twice  the  number  of  dominant  eigenvectors  if  the  tolerable  SNR 
degradation  is  to  be  no  more  than  3  dB.  The  numerical  experimi  nt  results  reported  by  Hung  et 
al.  [24]  agree  with  the  results  of  Kirsteins  and  Tufts  [29].  Even  if  PCI  is  used,  a  large  sonar  array 
may  still  have  the  number  of  sources  within  the  array  coverage  exceeding  the  number  of  available 
snapshots.  The  natural  question  is:  How  small  can  the  number  of  dominant  eigenvectors  be  as 
compared  to  the  number  of  sources? 

Weights  can  be  used  in  both  MUSIC  or  its  PCI  beamforming  equivalent  and  other  recent  sub¬ 
space  fitting  methods,  such  as  in  Owsley’s  Enhanced  Minimum  Variance  Beamforming  (EMVDR) 
[39,1,2].  (The  term  ’’weights”  used  here  refers  to  weighting  the  signal  eigenvectors.  It  is  not  the 
’’weight  vector”  used  in  the  better  known  MVDR  ABF  algorithm.) 

It  is  only  a  matter  of  time  before  someone  starts  to  apply  the  optimal  weights,  used  in  the 
WSF  for  the  parametric  signal  subspace  approach,  to  subspace  beamforming,  e.g.,  PCI  which  can 
be  considered  as  MUSIC’s  beamforming  variant. 

The  recent  papers  by  Farrier  and  Prosper  [16],  and  Farrier  [17]  are  very  important.  They 
managed  to  arrive  at  the  performance  of  those  jointly  estimating  DOA  algorithms  with  the  cost 
of  those  searching  one-at-a-time  spectrum  type.  They  made  a  more  efficient  modification  of 
Reiiiy,  Wong,  and  Reilly’s  Bayesian  approach  in  unknown  colored  noise  [44].  While  Reilly  et  al.’s 
approach  can  handle  100%-correlated  arrivals  and  colored  noise,  Farrier  and  Prosper  [16]  are  not 
far  off  as  one  of  their  methods  managed  to  resolve  two  98%-correlated  arrivals  in  one  scenario. 
By  using  this  same  method,  they  also  improved  the  noise  threshold  by  13  dB  in  another  case  as 
compared  to  MUSIC.  Among  the  three  weighting  schemes  evaluated  by  Farrier  and  Prosper  for 
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this  method,  the  optimum  weights  by  Viberg  and  Ottersten  outperformed  two  less  sophisticated 
ones  by  at  least  3  dB  in  noise  threshold.1 

The  word  beamformer  in  Farrier  and  Prosper’s  paper  is  a  misnomer  as  a  beamformer  produces 
a  time  series  at  each  grid  point  and  their  method  is  a  direction-finding  one  producing  a  null  spec¬ 
trum.  Consequently,  we  note  that  their  method  can  use  the  interpolation  scheme  by  Schmidt  [50]. 
There  he  shows  that  for  simple  scenarios  with  a  second  signal  at  least  a  Beamwidth  (BW)  away, 
the  grid  spacing  can  be  as  modest  as  0.25  BW  for  an  interpolation  accuracy  of  0.01  BW. 


2.3  Drawbacks  of  Oversimplifying  Stress  Measuring  Yardsticks 

It  turns  out  that  condition  number  is  not  the  only  stressing  factor  for  a  direction-finding  scenario 
and  there  are  quite  a  few  drawbacks  in  using  it  as  stress-measuring  yardstick.  Most  of  these 
deficiencies  are  manifested  under  the  following  special  cases  oitLci  singly  or  jointly:  the  two 
arrivals  are  100%  correlated,  of  equal  signal  strength,  and  with  temporal  phase  difference  0°  or 
180°.  Regardless  of  these  deficiencies,  it  is  still  a  useful  measure  summarizing  many  stressing 
factors  in  many  aspects. 

2.3.1  Infinite  Condition  Number  for  All  100%-Correlated  Arrival  Pairs 

When  the  two  arrivals  are  100%  correlated  so  that  the  smaller  eigenvalue  is  zero,  the  stress 
measured  by  condition  number  will  go  to  infinity.  Such  a  stress  measure  may  be  appropriate 
for  the  MUSIC  algorithm,  as  the  algorithm  fails  for  this  scenario,  but  not  appropriate  for  the 
scenario.  An  algorithm  like  weighted  subspace  fitting  can  handle  the  scenario  with  no  problem 
and  furthermore  asymptotically  realize  the  (finite)  stochastic  CRLB.  We  know  that  some  100%- 
correlated  scenarios  are  less  stressful  than  others,  and  in  these  cases,  the  efficiency  defined  by  the 
ratio  of  the  asymptotic  standard  deviation  of  the  parameter  estimation  to  CRLB  apparently  is  a 
better  stress  measure.  As  CRLB  is  available  for  all  scenarios,  one  would  naturally  wish  to  use  this 
as  a  universal  yardstick.  However,  there  is  a  practical  problem:  how  do  you  represent  all  possible 
scenarios  and  present  the  results  compactly?  So  far,  the  CRLB  is  presentable  only  for  uniform 
line  array,  and  even  then  not  compactly. 

2.3.2  Equal  Condition  Numbers  Corresponding  to  Different  CRLBs 

Apparently,  Schmidt  [49,  page  130]  was  the  first  one  to  use  the  noise- free  signal  eigenvalue  ratio  for 
stress  measurement  of  a  two-arrival  direction-finding  scenario.  He  recognized  that  the  ’’stress”  at 
the  combination  of  (£,  tj,  =  (0, 1, 1)  is  the  least  2  of  all  possible  ones.  The  two  signal-eigenvalues 
equal  the  normalized  unity  power  level  at  the  sensor  level.  We  believe  that  this  led  him  to  use 
an  expression  depending  only  on  the  eigenvalue  ratio  as  a  stress  measurement.  Evidently  he  was 
hampered  by  the  lack  of  an  effective  visual  aid  to  grasp  the  nontrivial  functional  dependence  of 
these  two  eigenvalues  on  all  possible  combination  of  the  product  of  the  two  unit-disks  representing 
the  temporal  and  spatial  correlation  between  the  two  arrivals.  This  possibly  led  him  to  conclude, 

‘These  less  sophisticated  weight  choices  include  the  one  that  Cadzow  used  in  his  recent  papers  [9.10]  where  he 
showed  that  his  eigenvalue  weighting  produced  results  marginally  better  than  the  ones  with  unit  weighting.  It  is 
somewhat  confusing  as  the  weight  choice  should  not  have  any  effect  when  there  is  only  one  dominant  eigenvector 
involved  such  as  in  the  case  of  two  100%-correlated  arrivals  used  in  these  papers.  If  the  dominant-eigenvectoi-only 
choice  was  abandoned,  the  rank-reduction  advantage  as  pointed  out  by  Ottersten  will  be  missed  too. 

2We  will  discuss  this  in  subsections  3.3.1,  3.3.3,  and  section  3.4. 
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without  proof,  a  statement  equivalent  to  one  that  the  eigenvalue  ratio  can  he  one  only  at  this 
(£,  Vi  ^ )  =  (0, 1, 1)  combination. 

We  will  see  that  when  the  two  arrivals  are  of  equal  strength,  i.e.,  ^  =  1,  the  two  eigenvalues 
are  equal  over  the  left  boundary  of  an  isosceles-right-triangle-like  region.  Furthermore,  the  two 
equally  dominant  eigenvalues  can  go  down  to  zero  jointly  so  that  the  condition  number  stays 
unity  throughout.  That  is,  the  condition  numbers  are  the  lowest  possible  there.  However,  both 
eigenvalues  tend  to  zero  as  r/  tends  to  zero,  i.e.,  the  baseline  of  the  region.  So,  when  the  noise  floor 
is  brought  in  for  example  in  the  symmetrical  multipath  scenario,  the  direction-finding  system  can 
still  be  stressed  significantly.  Yet  the  final  scenario  would  correspond  to  the  case  of  two  strong 
colocated  coherent  signals  exactly  equal  in  amplitude  but  totally  out  of  phase.  Even  though  the 
condition  number  remains  as  favorable  as  possible,  i.e.,  of  unity  value,  the  scenario  starts  from  the 
most  benign  one  in  stress  and  progressively  gets  more  and  more  stressful  and  eventually  gets  to 
the  point  of  being  totally  hopeless.  The  final  DF  problem  amounts  to  determine  the  value  of  two 
equal  numbers  when  it  is  only  known  that  the  difference  is  zero.  Because  the  two  signals  cancel 
out  each  other  exactly  at  the  colocated  emitters'  location  before  the  waves  can  propagate  out.  the 
sensor  measures  zero.  Therefore,  the  magnitude  of  largest  eigenvalue  in  unit  of  the  smaller  signal 
strength  at  the  sensor  level  is  a  good  complementary  stress  measurement. 

2.3.3  CRLB  Not  Appealing  to  Intuitive  Expectation 

The  effect  of  the  phase  difference  between  two  signals  on  direction-finding  is  apparently  not  a 
well-understood  subject.  Most  work  about  the  effect  of  phase  difference  in  direction-finding  is  for 
the  low-angle  radar  tracking  problem  with  two  closely  spaced,  coherent,  equipowered  plane-wave 
arrivals  impinging  on  a  UFA.  In  examining  candidates  for  the  stress  measure,  we  have  considered 
eigenvalue  ratio,  CRLB,  and  asymptotic  performance  results  of  some  subspace  algorithms  for  such 
scenarios.  However,  we  cannot  intuitively  reconcile  the  differences  among  these  candidates. 

If  the  center  of  the  array  is  chosen  as  coordinate  origin,  the  CRLB  is  symmetric,  is  maximum 
at  0°  or  180°,  and  is  minimum  at  90°.  The  analytical  derivation  of  this  result  is  available  in  the 
report  by  Evans,  Johnson,  and  Sun  [15].  When  we  have  two  100%-correlated  signals,  the  scenario 
points  in  the  isosceles  region  are  located  at  the  baseline  where  the  smaller  eigenvalue  is  zero. 
Furthermore,  when  the  two  signals  are  equipowered,  the  larger  eigenvalue  has  3  a  linear  variation 
varying  from  0  to  2  to  1  when  we  move  from  the  left-end  to  the  middle  and  then  the  right-end 
point  on  the  baseline.  We  try  to  use  CRLB  for  the  ULA  case  to  serve  as  the  stress  measure.  Note 
that  the  scenario  with  constructive  interference  has  the  same  CRLB  as  the  destructive  one,  i.e., 
it  is  symmetric  about  90°.  The  CRLB  has  its  minimum  occurring  at  this  90°  phase  difference 
point.  For  the  extreme  case  of  colocated  signals,  intuitively  one  would  think  that  the  constructive 
interference  case  is  a  more  favorable  scenario  than  the  destructive  one.  This  is  because  one  can 
decide  the  DOAs  from  the  constructively  combined  arrivals,  even  though  one  cannot  infer  whether 
more  than  one  arrival  contributes  to  the  measurement.  Neither  can  one  reconstruct  the  original 
two-component  amplitudes  given  oidy  the  final  total.  We  believe  this  dilemma  arises  because  of 
the  derivation  of  CRLB  from  a  rank-deficient  composite  steering  vector  matrix.  If  we  allow  the 
two  equipowered  arrivals  to  be  colocated  to  start  with  in  the  basic  data  model,  the  constructive 
scenario  may  be  differentiated  from  the  destructive  one  for  the  CRLB. 

For  the  equipowered  case,  the  destructive  interference  scenario  point  is  at  the  left  apex  of  the 
isosceles  region.  Both  eigenvalues  are  zero  and  they  are  global  properties.  I'he  eigenvalue  ratio 

at  this  point  is  of  the  0/0  form.  It  is  shown  in  subsections  3.1.1,  and  3.1.1  that  the  contours 

’We  will  discuss  this  in  subsection  3.3.2. 
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of  —  l  pass  through  this  point.  This  point  corresponds  to  the  total  cancellation  of  signals  at 
the  sensor  element  level.  The  total  cancellation  happens  not  only  at  the  center  of  the  array  but 
at  every  sensor.  The  two  steering  vectors  coincide  and  the  temporal  waveforms  negate  each  other 
completely.  Consider  two  strong  sources  behaving  this  way.  As  the  SNR  is  measured  by  turning 
off  one  source  at  a  time,  so  the  equal  SNRs  are  considered  as  high  and  CRLB  and  asymptotic 
results  apply. 

The  unified  framework  of  subspace-fitting  methods  leading  toward  Ottersten  and  Viberg’s 
WSF  formulation  for  this  two-source  problem  involves  different  weights.  The  optimum  weights 
for  the  high  SNR  case  involve  the  eigenvalue  ratio.  This  ratio  is  indeterminate,  and  the  asymptotic 
standard  deviation  can  assume  different  paths  leading  toward  the  anomaly  point. 

Tran  [58]  derived  and  plotted  the  2-d  mesh  surfaces  of  CRLB,  Det-ML,  WSF,  and  MUSIC 
against  phase  difference  and  SNR.  For  90%  temporal  correlation,  the  CRLB  surface  is  a  lower 
bound  for  the  other  asymptotic  surfaces  and  is  symmetric  with  respect  to  the  phase  difference  at 
90°.  Both  WSF  and  Det-ML  are  fairly  similar  to  each  other  because  the  SNR  used  is  high.  They 
lost  the  symmetry  at  0.2- BW  separation  but  retained  a  dip  with  180°  low-'  standard  deviation 
at  180°  than  at  0°.  MUSIC  at  0.2-BW  separation  is  monotonically  increasing  starting  from  180° 
passing  thru  90°  and  all  the  way  to  0°.  The  mesh  plots  show  that  the  symmetry  conditions  are 
present  at  1-BW  separation. 

We  need  to  reexamine  the  derivation  of  these  asymptotic  results  and  determine  whether  the 
rationale  in  dropping  high-order  terms  remains  valid  when  the  parameters  of  these  scenarios  are 
taken  into  consideration.  There  is  also  the  possibility  that  either  the  results  are  right  but  we 
need  some  interpretation  or  there  is  some  flaw  in  our  formulation.  Numerical  simulation  seems 
to  support  those  asymptotic  results.  Another  question  is  whether  high  SNR  is  relevant  under  the 
condition  of  perfect  cancellation. 

2.3.4  Redeeming  Value  for  Condition  Number 

Even  with  all  these  deficiencies,  condition  number  is  still  a  useful  measure  summarizing  stressing 
factors  in  many  cases.  When  the  phase  difference  between  two  signals  does  not  come  into  play  in 
the  eigenvalue  expressions,  such  as  when  the  two  arrivals  are  temporally  uncorrelated  or  spatially 
orthogonal,  eigenvalue  ratio  seems  to  be  a  plausible  stress  measure.  There,  a  large  eigenvalue  ratio 
will  reflect  the  effects  of  either  large  normalized  spatial-  or  temporal-  correlation  coefficient  or  high 
signal-strength  ratio  or  any  combination  of  them.  For  those  problems  with  100%-correlated  and 
closely  spaced  arrivals,  such  as  in  the  low-angle  (microwave)  radar  tracking  problem,  mean  square 
errors  or  CRLBs  may  be  substituted  for  condition  numbers  as  better  stress  measures,  and  the 
compact  presentation  for  scenario  parameter  tradeoff  need  be  used  with  cautions. 

Schmidt’s  original  MUSIC  algorithm4  uses  the  signal  eigenvalues  in  a  binary  manner.  That  is, 
a  hard  decision  between  noise  and  signal  is  made.  Viberg  and  Ottersten  avoided  making  a  binary 
decision  by  introducing  a  scale  equivalent  to  the  eigenvalue  ratios.  While  much  more  information 
may  be  gained  bv  using  eigenvalue  ratios,  we  are  still  constrained  by  the  information  that  is 
conveyed.  Other  measures  will  provide  us  more  infoi maliou  if  there  exists  a  manageable  way  to 

’Here,  we  include  the  MultiDimensional  (MD)-MUSIC  described  in  his  thesis  for  jointly  estimating  I)OAs 
simultaneously  which  was  motivated  by  considering  the  rank  deficiency  caused  by  100%-correlated  arrivals.  There, 
he  also  pointed  out  the  separable  nonlinear  least  squares  approach  for  MD-MUSIC.  The  name  MD-MUSIC  was 
coined  in  Roy’s  thesis  [d 5]  and  later  used  in  Viberg  and  Ottersten ’s  works.  While  we  have  occasionally  followed 
their  usage,  it  is  regretfully  a  misnomer.  Without  prior  exposure  of  its  definition,  people  tend  to  conjecture  that 
some  elevation  measurement  in  addition  to  that  of  azimuthal  one  is  involved.  A  better  replacement  is  "jointly 
estimating  MUSIC”  contrasting  the  original  one  of  searching  one-at-a-time  type. 
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present  it,  such  as  a  compact  map  showing  CRLB  or  asymptotic  performance  over  all  possible 
scenarios.  We  note  that  CRLB  depends  on  the  steering  vectors  not  just  through  the  function 
values  but  also  the  local  behavior  as  functions  of  the  direction-of-arrival  (DOA)  parameters.  This 
is  contrasted  to  the  real  triplet  (<f,  7,  ^  j  required  for  the  two  noise-free  eigenvalues  and  their  ratio. 
Therefore,  the  next  step  of  advancing  from  condition  numbers  to  CRLBs  is  a  very  significant  one. 

2.4  Special  Treatment  Given  for  Equal  Power  Arrivals 

The  following  are  the  reasons  for  giving  equal-power  arrivals  special  treatment: 

1.  It  happens  naturally  in  physical  phenomena,  e.g.,  for  the  case  of  low-angle  radar  tracking, 
see  Haykin  [22],  Reilly  and  Haykin  [41,42],  Haber  and  Zoltowski  [21],  and  Zoltowski  and 
Lee  [66].  The  problem  is  the  estimation  of  DOAs  of  a  radar  target  located  over  a  smooth 
reflecting  sea  when  the  main  beam  is  pointed  at  or  near  the  horizon. 

2.  Shown  in  subsection  3.4.2  is  that  for  given  ^  power  ratio,  the  A  slices  for  constant  phase- 
dependent  variable  £ s  are  parabolas.  The  equal-power  arrival  case  has  some  unique  analyti¬ 
cal  properties  stemming  from  the  fact  that  it  has  the  vertices  of  half  of  parabola  sections  fall 
within  the  range  of  parameters  of  interest.  Consequently,  there  are  infinite  slopes  involved 
here,  unlike  any  other  power  ratios  between  the  two  arrivals.  This  is  also  the  only  case  in 
which  it  is  possible  to  have  the  two  signals  to  cancel  each  other  completely. 

3.  Most  of  the  performance  characterization  using  either  Monte  Carlo  simulation  or  asymptotic 
analysis  include  equal-power  arrivals  as  examples.  For  instance,  see  Cadzow  et  al.  [6,7,8]. 
Kaveh  and  Barabell  [27],  Otterstcn,  Viberg,  et  al.  [34,36,37,59,60],  Pillai  and  R'won  [40,30], 
Reilly  and  Haykin  [41,42],  Stoica  and  Nehorai  [53,54,55],  Haber  and  Zoltowski  [21],  and 
Zoltowski  and  Lee  [66]. 

4.  The  concept  of  resolution  threshold  has  been  exploited  by  assuming  two  signals  of  equal 
power  impinge  on  a  ULA.  The  concept  was  started  by  Cox  [13]  for  CBF  and  ABF.  Kaveh 
and  Barabell  [27]  used  the  same  scenario  for  the  MUSIC  algorithm  for  incoherent  arrivals. 
Pillai  and  Kwon  [40]  extended  the  treatment  to  MUSIC  for  coherent  arrivals  via  spatial 
smoothing.  Kwon  and  Pillai  [30]  further  considered  the  resolution  threshold  issues  for  some 
special  cases5  of  three  equal-power  arrivals.  However,  the  concept  of  resolution  threshold 
has  outlived  its  usefulness  with  the  arrival  of  the  asymptotic  analysis  results  such  as  those 
from  Stoica  and  Nehorai  [53,54,55],  Ottersten,  Viberg,  et  al.  [34,36,59,60].  These  asymptotic 
analysis  results  yield  more  insight,  provide  much  more  information  enabling  one  to  make 
parametric  studies,  and  furthermore,  aro  applicable  to  general  sensor  array  configuration 
and  general  scenarios. 


'For  the  case  with  the  arrivals  equally  spaced  in  electric  angle,  the  required  threshold  SNRs  for  the  three-arrival 
case  to  resolve  the  three  are  about  twice  those  for  the  two-arrival  cases  in  dB  scale.  That  is,  the  presence  of  a  third 
source  which  is  symmetrically  located  with  respect  to  the  center  source  increases  the  SNR  required  to  resolve  the 
original  sources  approximately  by  a  factor  of  two  in  dBs. 

For  the  general  case,  Stoica  and  Nehorai  [55]  proved  that  the  asymptotic  Det-CRLB  increases  monotonically 
when  one  more  emitter  is  added  while  the  directions  of  and  the  signals  emitted  by  the  original  old  sources  remain 
the  same  Furthermore,  if  the  new  source  signal  is  uncorrelated  with  the  older  ones,  then  the  Sto-CRLB  increases 
monotonically.  They  stated  that  it  is  an  open  problem  whether  the  assumptions  of  uncorrelated  signals  might  be 
relaxed . 
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Chapter  3 


Eigenvalues  for  Two  Arrivals 


3.1  Eigenvalue  Equations 

We  follow  the  notations  in  Hudson’s  text[23,  pages  37-38,52-55]. 

(.)H  Hermitian  transpose 

(.)T  ordinary  transpose 

(.)*  complex  conjugation 

Re(.)  real  part  of  a  complex  number 

E(.)  expectation  operator 

the  quantity  for  the  new  coordinate 

j  V^T 

N  number  of  sensors  in  array 

x  iV-dimensional  array  data  vector,  x  =  Sm  +  n 

R  data  source  covariance  matrix,  E[xxH ] 

Si  steering  (column)  vector  for  the  ith  source 

5  composite  steering  vector  matrix,  S  —  {$,} 

<t>  generalized  cosine  between  and  S2 

p  normalized  correlation  between  first  and  second  baseband  signals 

tti  signal  power  for  the  stronger  source 

7T2  signal  power  for  the  weaker  source 

Pi  ith  eigenvalue  of  R 

Ai  larger  eigenvalue  of  the  noise-free  part  of  R,  normalized  to  Nn 2 

A2  smaller  eigenvalue  of  the  noise-free  part  of  R,  normalized  to  Nx? 

tj  (1  —  jd>|2)(l  -  |p|2),  independent  of  the  angle  difference  between  the  normalized 

temporal  and  spatial  correlation  coefficients  p  and  4>,  both  constrained  to  unit  disks. 

£  Re(p4>*)  =  \p\\<f>\  cos(argp  -  arg<£),  dependent  on  angle  difference 

a  the  coefficient  vector  expressing  an  eigenvector  as  a  linear  combination  of  steering  vectors 

z,  position  vector  of  the  ith  sensor  element 

k  unit  direction  cosine  vector  for  the  direction  of  arrival(DOA)  of  a  signal 

c  speed  of  wave  propagation 

6  the  difference  between  two  coordinate  origins 

u>o  the  center  frequency  for  the  narrowband  signals 
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(xn,yn)  xn  =  (n  -  1  )d,  yn  =  0,  the  elements’  location  of  Uniform  Line  Array  (ULA) 
d  the  uniform  interelement  spacing 

a,ao  directions  of  two  arrivals  relative  to  broadside 

6  sin  ct  —  sinao,  the  reduced  angular  variable 

k  2?r /wavelength,  wave  number 

kd$  the  electrical  phase  angle  from  element  to  element  along  the  ULA 

4.dB-down.bwJnj9  4  dB  down  beamwidth 

0 ]T<icJ>w  separation  between  two  arrivals  in  fractional  beamwidth 

In  complex  notation,  the  N -dimensional  sensor  output  vector  x  is  modeled  by  a  linear  combi¬ 
nation  of  steering  vectors  s ;  of  M  sources  plus  an  IV-dimensional  noise  vector  n 

x  =  Sm  +  n, 

where  mi  is  the  complex  amplitude  of  the  /-th  arrival.  The  usual  assumptions  of  white  noise  and 
zero  correlation  between  each  signal  and  noise  lead  to 

R  =  E[xxh]  =  S  E[mmH ]  SH  +  cr2/, 

where  a 2  is  the  variance  of  the  spacial  white  noise.  The  eigenvalue  p  and  eigenvector  u  of  R  is 
defined  by 

Ru  =  hu- 

Hudson  in  [23]  used  the  following  developments  restricted  to  the  two-arrival  case.  The  eigen¬ 
values  m  of  R  are  related  to  the  eigenvalues  Aj  of  E[mmH]  SH S  through  the  relation 

Hi  =  K  +  o2. 

The  general  case  was  mentioned  in  Kwon  and  Pillai[30].  These  known  results  are  rederived  in  the 
following.  We  start  from 

(5  E[mmH]  SH  +  a2I)u  =  /ru, 
or 

( S  E[mmH ]  SH)u  =  (h  -  o2)u  =  Xu. 

Therefore,  the  non-noise  eigenvectors  are  linear  combinations  of  the  steering  vectors  of  the  DOAs, 

u  =  5(  %.»1  S«u) ;  Sa 

A 

The  non-noise  eigenvalue  can  be  derived  from 

S(  E[mmH]  SHS  -  X I)a  =  0. 

With  the  usual  assumptions  that  the  steering  vectors  for  different  DOAs  are  linearly  independent, 
we  have 

(  £[mm^]  SH S  -  XI)a  =  0. 

Note  here,  the  components  of  each  eigenvector  in  this  new  non-Hermitian  eigenproblem  correspond 
to  the  coefficients  expressing  the  corresponding  eigenvector  of  the  original  Hermitian  eigenproblem 
expressed  as  a  linear  combination  of  steering  vectors. 
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In  general,  the  product  of  two  Hermitian  matrices  is  not  Hermitian.  In  fact,  for  the  two-arrival 
case,  we  can  see  from  the  product  matrix  in  Hudson’s  text[23,  pages  52-55]  that  the  diagonals  are 
no  longer  real  and  the  off-diagonal  entries  are  not  complex  conjugate  of  each  other.  However,  the 
eigenvalues  of  this  product  matrix  should  be  real  and  non-negative  as  they  are  from  the  original 
Hermitian  eigenproblem.  Therefore,  the  sum  of  the  diagonal  entries  are  still  real  because  of  the 
invariant  property  of  the  trace  after  orthogonal  transformations.  By  using  a  similar  but  simpler 
factoring  argument,  Speiser  and  Arnold[51]  in  fact  derived  useful  bounds  for  singular  values  of 
direction-finding  problems  for  the  general  case.  As  the  emphasis  of  this  report  is  on  two  arrivals, 
we  will  evaluate  his  bounds  in  section  3.6  in  this  context. 

We  return  to  the  consideration  of  two- arrival  noise- free  dominant  eigenvalues.  These  two 
eigenvalues  for  the  eigenproblem 


where 

*?  =  (!-  I^X1  -  M2)> 

£  =  R e(p<t>’)  =  \p\\(t>\  cos(argp  —  arg(£). 
Let  7Ti  >  tt2.  The  eigenvalue  problem  can  be  rewritten  as 


Henceforth,  the  eigenvalues  will  be  redefined  as  their  normalized  values  with  respect  to  Nv-i,  so 
that 


Before  finishing  the  introduction  of  symbols  and  definitions,  we  point  out  the  importance  of  the 
sensor  coordinate  origin.  When  two  arrivals  are  100%  correlated,  the  phase  difference  between  the 
two  source  signals  plays  an  important  role.  The  specification  of  this  variable  depends  on  the  choice 
of  the  coordinate  origin.  Frequently,  performance  such  as  CRLB  is  plotted  against  this  variable. 
If  two  different  origins  are  used,  the  same  results  can  be  seen  differently  but  are  equivalent  with  a 
horizontal  shift  and  possibly  with  the  help  of  27r-modulo  wrap-around.  Naturally,  the  covariance 
matrix  and  their  eigenvalues  are  not  dependent  on  the  choice  of  the  origin,  but  the  expressions 
for  the  phase  angles  of  both  the  steering  vectors  and  complex  source  amplitudes  of  the  signals 
are.  The  following  make  these  discussions  clearer. 

Denote  the  difference  between  the  two  coordinate  origins  by  6,  i.e.,  the  sensor  position  de¬ 
scriptions  in  these  two  coordinate  systems  are  related  by  z*  =  z,  +  6  for  all  sensor  elements.  The 
signal  part  of  the  sensor  data  vector  at  the  ith  sensor  is 

2  2 

(5m),  =  £  Slkmk  =  £  ej^^cmk 

*=l  k=l 
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2  2 

=  £  (e^ZVc)  (eiu *«ls/cmk>i  =  £  §tkmk 

k=l  k=\ 

=  (5m),. 

Define  D  =  diag{e^K^s'lcmk).  Then  we  have  5  =  SD~1,  and  m  =  Dm.  Note 
(  E[mrhH] )«  =  (D  £(mm"]  DH)«  =  ew»(«.-*i)T^  (  £[mm"]  )t/, 

and 

(5wS)i,  =  ( D~HSHSD~l)u  =  (5w5)i,. 

The  angular  difference  between  the  temporal  and  spacial  correlation  coeffients  is  therefore  inde¬ 
pendent  of  the  coordinate  choices,  as 

arg(  E[rhmH]  ),/  -  arg  (SHS),i  =  arg(  E[mmH]  ),,  -  arg(5H5),7. 

The  assignment  of  an  array  coordinate  origin  synchronizes  the  angles  of  zero  references  of  the  two 
disks. 

In  the  special  case  of  plane  waves,  an  array  is  called  pairwise  symmetric  or  inverse  symmet¬ 
ric  [13]  if  for  each  sensor  located  away  from  the  origin  at  position  coordinates  (x,,  y,-,  Z{)  there 
is  also  a  sensor  located  at  —  j/,-,  —  z,-).  The  inner  product  of  two  steering  vectors  for  plane 

waves  impinging  onto  a  pairwise  symmetric  array  is  real.  For  such  case,  the  usual  two-dimensional 
unit-disk  associated  with  spatial-correlation  <f>  degenerates  into  a  one-dimensional  one.  Then  the 
angle  difference  between  the  temporal-  and  spatial-correlation  coeffients  is  reduced  to  the  temporal 
phase  difference  between  the  two  arrivals  measured  at  the  origin  of  the  symmetric  array. 

Even  though  such  relation  is  derived  from  plane  waves  impinging  onto  pairwise  symmetric 
arrays  for  convenience,  we  call  r)  the  phase-independent  variable  and  £  the  phase-dependent  vari¬ 
able. 
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3.2  An  Isosceles-Right-Triangle-Like  Region 

One  can  see  easily  that  the  normalized  As  and  their  ratio  A1/A2  are  functions  of  the  three  param¬ 
eters  77,  £,  and  X1/X2.  It  is  clear  that  0  <  77  <  1.  But  77  and  £  are  not  entirely  independent  of 
each  other.  In  the  following,  we  will  show  that  given  77,  £  is  bounded  depending  upon  the  77  value 
which  in  turn  depends  only  on  \p\  and  |<£|. 

As  £  is  the  triple  product  of  |p|,  |<£|,  and  cos(arg/>  —  argd>),  it  suffices  to  bound 

|£|max,  given  r ;  =  IpII^I 

subject  to  the  equality  constraint 

(i  -  M!)(i  -  W2)  =  v. 


Using 


( \pM )2  =  =  1 


1 

=  -(l-H2)  +  (1  +  7,)- 


1  -  \p? 


1  -  \p\2' 


and  setting  to  zero  the  derivative  of  the  above  expression  with  respect  to  (1  —  |p|2)  to  maximize 
£,  we  have 

—  1  + - - - =  0 

fl  -  lol2  V  ’ 

V*  lrlm«x  (,  given  r)J 

then  we  have 

IHmax  given  *7  ~ ^  ~  \/Vw 

From  the  definition  of  77,  a  similar  relation  holds  for  <f>,  i.e., 


I^lmax  4,  given  v  ^  y/V- 

Therefore,  we  have 

(IHI^Dmax  given  rj  =  1  ~  y/V- 

For  given  n\/n 2,  the  domain  in  the  77  and  £  plane  over  which  As  and  their  ratio  Aj/A2  are 
defined  is  the  isosceles-right-triangle-like  region  bounded  by  one  line  and  two  parabolas: 

77  =  0,  and  77  =  (£  ±  l)2. 


We  note  the  zero  slope  of  each  parabola  where  it  meets  the  straight  line  and  the  90°  angle  that 
the  two  parabolas  intersect  each  other. 

We  will  see  later  that  among  the  three  parameters  (£,77, X1/X2)  appearing  in  the  eigenvalues’ 
expressions,  77  is  the  most  important  one.  In  the  following,  we  use  figure  3.1  to  provide  some 
characterization  of  the  many-to-one  mapping  from  the  Cartesian  product  of  the  two  unit-disks 
onto  the  isosceles-like  region  with  respect  to  77.  We  note  that  for  given  77,  max(|p|,  |</>|)  =  —  y/rj. 

Therefore,  the  apex  point  (0,1)  corresponding  to  77  =  1  can  only  come  from  the  two  centers  of 
the  unit-disks,  i.e.,  p  —  <j>  =  0.  As  77  decreases  from  1  to  0,  i.e.,  we  move  from  the  apex  towards 
the  baseline,  we  shade  the  part  of  the  Cartesian  product  of  the  two  unit-disks  which  is  mapped 
into  the  correspondingly  shaded  area  of  the  isosceles-like  region.  Therefore,  as  the  shaded  area  of 
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the  isosceles-like  region  expands  from  the  apex  towards  the  baseline,  the  shaded  areas  of  the  two 
unit-disks  expand  from  the  centers  toward  the  unit  circumferences. 

When  we  move  one-fourth  the  way  toward  the  baseline,  i.e.,  77  >  0.75,  we  find  max(jpj,  |<£j)  = 

\J  1  —  \/0-75  ss  0.36.  Let  us  shade  the  corresponding  range  and  domain  of  this  many-to-one 
mapping  using  45°  hatching  lines  for  this  i)  >  0.75  region.  Similarly,  when  we  have  moved  halfway 
toward  the  baseline  in  the  isosceles-like  region  t 7  >  0.50,  we  find  max(|/>|,  |<£|)  =  \J  1  -  \ATE  « 
\/0.3  as  0.54.  Let  us  overlay  onto  the  earlier  picture  by  shading  the  corresponding  range  and 
domain  using  offset  45°  hatching  lines  for  this  rj  >  0.50  value.  We  note  that  the  previously  45° 
hatched  region  are  double-hatched  now.  For  any  point  in  the  incremented  area  of  the  isosceles-like 
region  corresponding  to  the  singly  45°  hatched  but  not  double-hatched  area,  i.e.,  0.75  >  77  >  0.50, 
at  least  one  of  the  originating  p  or  (f>  must  be  located  in  the  correspondingly  shaded  incremental 
annulus,  i.e.,  0.36  <  |p|  or  |<£|  <  0.54. 

We  also  note  that  77  =  (1  —  |<£|2)(1  -  jp|2)  is  the  height  of  the  unshaded  trapezoidal-like  region. 
We  can  associate  X77  as  the  product  of  the  areas  of  the  two  unshaded  annulus  regions,  x(l  —  |(/>|2) 
and  x(l  —  |/>|2).  As  this  height  reduces  from  1  to  0,  i.e.,  the  unshaded  trapezoid  decreases  from 
the  original  whole  isosceles  all  the  way  to  the  zero  thickness  baseline,  the  corresponding  white 
annulus  shrinks  from  the  original  unit-disk  to  a  zero  thickness  circle. 

When  the  scenario  is  stressful,  r\  is  small,  the  corresponding  annuli  with  unit  outer  radius  are 
thin.  The  area  of  the  temporal  correlation  p  annulus  can  be  approximated  by 

*(i  -  M2)  =  *(i  +  M)(i  -  M)  *  2*-(i  -  |H). 

Suppose  the  temporal  correlation  coefficient  has  a  magnitude  of  99%  and  up,  then  the  thickness 
of  the  p  annulus  is  0.01,  the  area  of  the  annulus  can  be  approximated  by  27t(0.01)  by  letting  the 
mean  circumference  of  the  annuls  assuming  the  outer  circumference  of  the  unit-disk,  2 n.  The 
approximation  is  99.5%  accurate.  If  the  spatial  correlation  coefficient  has  a  magnitude  of  90% 
and  up,  then  the  thickness  of  the  </>  annulus  is  0.1,  the  area  of  the  annulus  can  be  approximated 
by  2x(0.1).  The  approximation  accuracy  is  95%.  We  see  that  a  combination  of  the  high  temporal 
and  spatial  correlations  yield  an  77  of  approximately  (—20  +  3)  +  (  — 10  +  3)  =  —24  dB. 

For  given  77,  we  are  interested  in  the  maxima  of  |p|,  |<£|,  and  |£|.  Given  77,  77  =  (1  -  |<£|2)(1  - 
|p|2)  <  min(l-|p|2,  l-|(^>j2).  That  is,  max(|p|2, |$|2)  <  1  —  77,  or  max(|p|,  |<^>|)  <  y/1  -  rj.  Therefore, 
max(|7>|)  =  y/i  -  77  occurs  when  1  -  \<f>\ 2  =  1  or  <f>  =  0.  Similarly,  max(|<£|)  =  ^1-77  occurs  when 
p  =  0.  Both  scenario  points  are  at  the  vertical  axis,  £  =  0. 

Consider  77  =  |,  i.e.,  three-fourth  of  the  way  toward  the  baseline.  |£|mjn  given  77  =  ^  ocn,rs 

at  the  vertical  axis,  with  |p|  or  \<t>\  =  i/l  -  77  =  ^  %  0.866.  We  already  know  that,  from  the 
derivation  of  the  isosceles  region,  |£|max  given  77  occurs  at  IH  —  l<£|  =  -  \/V  —  j  ~  0.707, 

i.e.,  the  two  parabolic  boundaries.  These  are  summarized  in  table  3.1. 

Figure  3.1  need  to  be  interpreted  with  cautions.  The  domain  of  the  many-to-one  mapping, 
a  subset  of  the  4-d  Cartesian  product  of  two  disks,  for  a  range  having  the  shape  of  a  shaded 
2-d  subtriangular  region  cannot  be  explained  briefly.  The  labels  and  sizes  of  the  shaded  disks  in 
figure  3.1  need  be  interpreted  according  to  table  3.2. 

Next,  we  use  figure  3.2  to  characterize  some  special  positions  in  the  isosceles  region.  Hudson’s 
textbook  presented  4  out  of  the  17  characterizations  shown  there.  These  four  cases  are  denoted 
by  the  asterisks. 
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Table  3.1:  Expressions  of  |£|,  |p|,  and  |<£|,  given  t)  = 


ICI 

w 

min  =  0  (vertical  axis) 

0  (center  of  disk) 

max=  y/l  -  T)  =  ^  a  0.866 

max  =  y/\  —  tj  =  ^  «  0.866 

0  (center  of  disk) 

max  =  1  -  y/rj  =  |  =  0.500 
(parabolic  boundary) 

#«0.707 

yjl-y/ij  =  #  «  0.707 

Table  3.2:  Interpretations  of  the  <f>  unit-disk  for  two  uncorrelated  plane  wave  arrivals  impinging 
onto  a  uniform  linear  array. 


V 

max|<£|2  =  l-r? 

max|<£|  =  y/l  -  Tj 

Separation  Between  Arrivals 

1 

0 

0 

nulls,  i.e.,  1,2,-  •  -,-BW 

3 

4 

l 

4 

\  =  0.500 

>i-BW 

1 

2 

1 

2 

i-BW 

1 

4 

3 

4 

1K11I1M 

<  i-BW 

0 

1 

i 

coincident  mainlobes  or  grating  lobes 
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3.3  Special  Cases  at  the  Apex,  the  Baseline,  and  the  Vertical 
Axis 

We  now  characterize  some  special  cases  of  interest.  They  are 

1.  at  the  apex,  i.e.,  (£,77)  =  (0,1), 

2.  at  the  baseline,  i.e.,  77  =  0, 

3.  at  the  vertical  axis,  i.e.,  £  =  0. 

3.3.1  Special  Case  at  the  Apex,  i.e.,  (£,77)  =  (0,1) 

Here  we  have  the  two  arrivals  0%  correlated  and  the  two  steering  vectors  orthogonal  to  each  other 
simultaneously.  Note 

X\  —  — ,  and  X2  -  1  or  101og10(A2)  =  0  dB. 

X2 

This  is  the  most  benign  scenario  for  a  direction-finding  system  where  the  two  arrivals  behave 
essentially  independently  of  each  other. 


3.3.2  Special  Case  at  the  Baseline,  i.e.,  77  =  0 

This  is  the  case  where  the  source  amplitudes  of  the  two  arrivals  are  100%  correlated  or/and  the 
two  corresponding  signal  vectors  are  complex  scalar  multiples  of  each  other  or  both. 

We  have  the  following  results: 


A2  =  0,  or  101og10(A2)  =  -00  dB, 


TTl 


Ai  =  —  +  1  +  2{. 


7T2 


We  note  that  when  £  changes  from  -1  to  0  to  +1,  Aj  varies  linearly  in  £  with  the  slope  equal  to 


V* from 


(\/~  “  1)  to(-i  +  l)to(./-i  +  l) 


7T  2  7T2 

Therefore,  the  two  vertices  corresponding  to 


(<,»?)  =  (-1,0)  and  (1,0) 

are  associated  with  totally  destructive  and  constructive  interferences  respectively.  Here  the  two 
arrivals  are  colocated  with  respect  to  the  steering  vectors,  e.g.,  along  the  same  direction  of  arrival, 
and  the  source  amplitudes  are  exactly  opposite  in  sign  or  equal  to  each  other.  At  the  mid-point, 
(tiV)  =  (0,0),  we  have 

Ai  =  —  +  1. 

7T2 

When  the  two  100%  correlated  arrivals  are  ’’quadraturely  related,”  this  is  the  place  where  it  occurs. 
(In  direction-finding,  any  two  arrivals  are  called  coherent  as  long  as  their  source  amplitudes  are 
100%  correlated,  even  if  their  phase  difference  is  90°. ) 
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3.3.3  Special  Case  at  the  Vertical  Axis,  i.e.,  £  =  0 
We  have  any  of  the  following  three  possibilities: 

1.  p  =  0,  i.e.,  the  source  amplitudes  of  the  two  arrivals  are  0%  correlated, 

2.  4>  =  0,  i.e.,  the  two  steering  vectors  are  orthogonal  to  each  other,  e.g.,  a  delay  and  sum  beam 
for  the  case  of  ULA  in  the  direction  of  the  first  arrival  will  produce  a  null  in  the  direction 
of  the  second, 

3.  cos(argp  —  arg<£)  =  0,  i.e.,  the  2  phasors  on  the  two  unit-disks  are  perpendicular  to  each 
other. 


We  have 


(ff  +  1)2 


This  can  be  rewritten  in  a  more  revealing  form  of  a  parabola  in  A  and  rj 


5<?  +  i) 

Z  7T2 


_  *1  _  _*i 

7T2  7T2 


V- 


(£  +  l)2 


4^ 

*2 


This  equation  describes  that  the  two  curves  of  Ai  and  A2  as  functions  of  77  for  £  =  0  pertain  to 
two  parts  of  the  same  truncated  parabola  in  A  and  77.  Because  the  arithmetic  mean  of  ^  and  1  is 
never  smaller  than  the  geometric  mean,  the  vertex  of  this  single  parabola  lies  outside  the  domain 
of  interest  except  when  ^  =  1.  In  the  exceptional  case,  there  is  no  truncation  of  the  parabola  in  A 
and  77  and  the  two  curves  meet  at  unit  height.  As  77  decreases  from  1  to  0,  the  concave  downward 
curve  of  Ai  increases  from  to  ^  -f  1,  whereas  the  concave  upward  curve  of  A2  decreases  from 
1  to  0. 

Taking  partial  derivatives  of  the  last  equation  of  the  two  curves  with  respect  to  77  and  rear¬ 
ranging,  we  have  for  7  =  1,2 


-ILL 

”2 


At  the  baseline  where  77  =  0,  Ai  =  +  1,  and  A2  =  0,  we  have 


£ai(o,o,  =  -£*o.o 


=  -1  + 


(£  +  1) 


At  the  apex  where  77  =  1,  Aj  =  and  A2  =  1,  we  have 


’Ll 

”2 


=  -1  - 


(ft +  *)' 


For  large  and  as  77  ranges  from  1  to  0,  the  theoretical  decrease  in  magnitude  of  the  slope 
of  the  slightly  concave  downward  curve  of  Aj  is  hardly  noticeable.  For  all  practical  purposes,  the 
A2  curve  is  a  straight  line  with  slope  1  and  values  between  1  and  0,  and  that  of  Ai  is  a  straight 
line  with  slope  -1  between  1  and  2.  When  scaling  is  used,  the  A!  and  A2  curves  plotted  together 
will  be  fairly  flat  for  large  *l  values. 
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When  ^  =  1  and  as  77  ranges  from  1  to  0,  the  slope  of  the  concave  downward  curve  of  Ai 
increases  from  — oo  to  -1/2,  whereas  the  slope  of  the  concave  upward  curve  of  A2  decreases  from 
4*oo  to  1/2.  The  ±oo  values  of  the  slopes  of  the  A2  and  Ai  curves  with  resper.,  to  7  at  77  =  1 
indicate  the  rapid  rate  of  the  increase  in  ’’stress”  to  the  direction-finder  for  the  case  of  ^  =  1. 
However,  we  also  note  that  the  "stress”  at  this  combination  of  (£,77,  ;£)  =  (0, 1, 1)  is  the  least  of  all 
possible  ones.  Unless  special  care  is  given,  commonly  used  plots  will  not  reflect  the  infinite-slope 
case  well. 

As  we  try  to  use  the  eigenvalue  ratio  ^  as  a  stress  indicator,  the  rest  of  this  subsection 
concerns  the  ratio  ^  as  a  function  of  and  77.  We  separately  discuss  the  two  cases  of  jp-  >  1 

and  «  1. 

*2 

When  >  1,  we  first  approximate  the  radical  expression  as 


1 


1  - 


(£  +  l)2 


Then  the  ratio  ^  is  given  approximately  by 


2r? 

2L 

*2 


'2 


That  means  the  stress  from  7?  and  ^  is  multiplicative  (divisive)  in  the  direct  scale.  In  dB  scale, 
this  is  additive,  i.e., 


When  ^  «  1,  we  have 

IT'i 


^(in  dB)  %  (  — )(in  dB)  -  i/(in  dB). 
A  2  7T2 


A^  _  1  4-  y/1  -  77  _  (1  4-  v/1  -  T])Z  _  2  -  7?  4-  2y/  1  -  77  ^  4 
A2  1  -  \/\  —  77  1  -  (1  -  77)  77  ""77' 

We  used  tj  <C  1  in  the  last  step.  Translated  into  dB  scale,  this  is 

-r^(in  dB)  ss  —  7/(in  dB)  4*  6. 

A2 

This  extra  6  dB  accounts  for  the  presence  of  the  infinite  slope  of  A  with  respect  to  tj  at  the  apex 
for  the  equal-power  arrival  case. 

These  analytical  results  will  be  verified  in  figure  5.29  and  discussed  in  subsection  5.4.1. 


3.4  Special  Case:  Equal  Strength  Arrivals,  i.e.,  tti/ tt2  =  1 

The  normalized  eigenvalues  assume  the  following  expressions 

*»  =  (!+£)  l±  \jl  ~  li  =  (1  +  \/(!  +  02  ~  V- 

We  note  that 

A, (0,1)  =  A2(0,1)  =  1,  and  A,(-1,0)  =  A2(£,0)  =  0. 

The  first  equation  states  that  the  two  surfaces  Aj(£,77)  and  A2(^,  77)  coincide  at  the  apex  of  the 
isosceles-right-triangle-like  region,  as  they  should.  The  second  equation  states  that  both  eigen¬ 
values  coincide  again  and  vanish  at  the  left  vertex  of  the  baseline  corresponding  to  destructive 
interference. 

The  two  eigenvalues  are  equal  (not  necessarily  zero)  only  if  what  is  inside  the  radical  vanishes. 
For  the  case  Tj/jr2  =  1,  the  resulting  equation  coincides  with  that  of  the  left  parabola  in  the  £ 
and  77  plane  bounding  the  region  over  which  the  As  are  defined.  The  common  value  they  share  is 

A,(£,(1  +  02)  =  (1  +  0  =  V^ 

However,  we  caution  that  unless  special  effort  is  made  in  plotting  results,  the  two  displayed 
surfaces  or  curves  may  not  always  meet  at  the  supposed  places.  This  is  because  the  boundary 
curves  defining  parabolas  in  the  £  and  77  plane,  in  general,  do  not  pass  through  the  grid  points 
used  in  the  (  and  77  plane. 

At  the  baseline  77  =  0,  we  have  A2  =  2(1  +  £).  This  means  that  when  £  changes  from  -1  to  0 
to  +1,  A!  varies  linearly  in  £  with  unit  slope  from  0  to  2  to  4.  At  the  vertical  axis  £  =  0,  we  have 
the  following  parabola  forms 

Aj  =  1  ±  \J\  —  77. 

3.4.1  Straight  Line  Contours  for  Constant  Eigenvalue  As 

We  now  investigate  the  contour  plots  of  these  two  surfaces  by  using  the  following  rearranged 
equations 

a,-(i  +  o  =  ±v/(i  +  <02-t?, 

or 

77  +  A  2  =  2Ai(l  +  £). 

These  contours  are  straight  lines  in  the  £  and  77  plane  with  positive  slopes  2A,  and  intercept  the 
vertical  line  passing  through  the  left  vertex  (-1,0)  of  -  A?.  The  mesh  and  contour  plots  of  and 
A2  are  shown  in  figures  3.3  and  3.4  respectively. 

As  the  contour  values  of  Ai  increase  from  0  to  2  to  4,  the  slope  increases  linearly  from  0  to  4 
to  8,  and  the  intercept  decreases  quadratically  from  0  to  -4  to  -16.  These  contours  intersect  the 
baseline  in  equal  intervals  as  we  noted  early. 

As  the  contour  value  of  A2  increases  from  0  to  1,  the  slope  increases  linearly  from  0  to  2,  and 
the  intercept  decreases  quadratically  from  0  to  -1.  Note  the  0  contour  coincides  with  the  baseline. 


31 


3.4.2  Parabolic  Slices  for  Constant  Phase-Dependent  Variable  £s 


For  constant  £  slices,  the  characteristic  equation  can  be  rearranged  to  that  of  a  parabola  in  A  and 
77  as 

[A-(i+£)j2  =  -[r?-(i  +  £)2  • 

Figures  3.5  plots  for  Ai  And  A2  vs  77  At  £=-0.50,-0.25,  0,  +0.25, +0.50  For  =  0  dB. 

For  constant  £  <  0  slices,  the  vertexes  of  the  vertica1  parabolas  in  A  and  77  will  be  exactly 
over  the  parabola  77  =  (£  +  l)2  in  the  horizontal  £  and  77  plane  defining  the  left  boundary  of  the 
isosceles-right-triangle-like  region,  with  heights  of  ±(1  +  £)  =  ±y/rj  as  described  earlier  from  a 
slightly  different  perspective.  The  ±00  values  of  the  slopes  of  the  A2  and  Ai  curves  with  respect 
to  77  at  the  vertexes  of  these  parabolas  indicate  the  rapid  rate  of  the  increase  in  "stress”  to  the 
direction-finder  for  the  case  of  ^  =  1.  For  constant  £  >  0  slices,  the  vertical  parabolas  in  A  and  77 
will  have  their  vertexes  located  over  the  continuation  of  the  just  mentioned  left  boundary  parabola 
in  the  horizontal  £  and  77  plane.  Therefore,  the  apexes  are  outside  the  isosceles-right-triangle-like 
region,  and  furthermore,  beyond  77  >  1.  The  curvature  of  the  parabola  is  1  at  the  vertex  and  is 
independent  of  the  £  value;  i.e.,  all  the  parabolas  are  translational  copies  of  each  other. 

Taking  partial  derivatives  of  the  parabolic  equation  with  respect  to  77  and  rearranging,  we  have 


d  -1  ±1 

dr,  2  [A  —  (£  +  1)]  “  \/(l  +  £)2  —  77' 


At  the  77  =  0  baseline,  we  have 

For  the  following  three  special  points,  they  are 


±1 

2(1  +£)' 


d_ 

dr. 


A(-1,0)  =  ±00, 


|ia,-.,o)=4 


We  note  that  the  ±00  values  are  not  very  useful  mathematical  expressions,  as  only  a  single  point 
is  defined  in  the  77  direction  at  the  left  vertex  as  the  domain  of  interest  is  concerned.  Similar 
remarks  hold  for  the  right  vertex.  These  values  suggest  the  neighborhood  behavior  about  the 
vertexes,  however. 

For  given  £  and  ^  =  1,  as  for  all  other  values,  the  stress  of  a  direction-finding  scenario  is 
highest  when  the  baseline  77  is  approached.  We  are  interested  in  the  local  behavior  at  the  stressful 
region  about  77  =  0.  In  this  region,  the  large  eigenvalue  Aj  will  not  change  much  in  terms  of  its 
order  of  magnitude.  This  is  not  the  case  for  the  small  eigenvalue  A2,  as  it  is  the  difference  of  two 
nearly  equal  numbers.  The  local  first-order  approximations  to  Ai  and  A2  near  77  =  0  are  therefore 


M£,I7)«  A,(£,0)  =  2(1  +£), 


A2(£,  77) 


1 

2(1  +  £)' 


Near  77  =  0  and  at  the  three  particular  £  values  of  -1,  0,  and  1,  we  have 


A2(  —  1 , 77)  «  0077,  A2(0, 77)  =s  A2(l,T7)ss|. 

We  can  confirm  these  analytical  findings  by  reviewing  figure  3.5.  These  discussions  about  the 
local  behavior  of  the  eigenvalues  near  the  baseline  also  serve  as  a  preview  of  the  next  figure  3.6 
in  the  next  subsection. 
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lambda_l  &  2  t,  lambda  1  &  2 


igure  3.5:  Parabolic  slices.'Aj  and  A2  vs  r\  at  £=-0.50,-0.25,  0,  +0.25, +0.50  for  ^  =  0  dB. 


Figure  3.6:  Hyperbolic  slices:Ai  and  A2  vs  £  at  rj  =  0.0,0.25,0.50,0.75,1.00  for  ^  =  0  dB. 
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3.4.3  Hyperbolic  Slices  for  Constant  Phase-Independent  Variable  77s 

For  constant  77  slices,  the  characteristic  equation  can  be  rearranged  to  that  of  a  hynerbola  in  A 
and  £  as 

A2  -  2£A  -  2A  +  77  =  0. 

Figures  3.6  plots  for  Ai  and  A2  vs  77  at  77  =  0.0,0.25,0.50,0.75,1.00  for  ^  =  0  dB. 

We  defer  the  derivation  and  characterization  of  1  ese  hyperbolas  to  appendix  A  and  the  next 
section  for  the  general  case.  We  mention  a  few  results  here. 

The  two  straight  lines  represented  by  Ai  and  A2  over  the  baseline  77  =  0  are  those  of  a 
degenerate  hyperbola  from  a  one-parameter  77  family  of  hyperbolas  in  A  and  £  sharing  common 
asymptotes  when  the  parameter  77  is  0.  The  common  asymptotes  are  the  two  straight  lines. 

Having  defined  the  common  asymptotes  of  this  family  of  hyperbolas,  the  next  important 
property  is  how  the  vertex  varies  when  77  increases  from  0  to  1.  Define  Ey  as  the  distance  from 
the  vertex  to  the  intersection  of  the  common  asymptotes.  We  note  that  Ey  is  proportional  to  y/rj, 
i.e.,  Ey  at  77  =  |  is  half  of  Ey  at  77  =  1. 


3.4.4  Contours  for  Eigenvalue  Ratio 

Consider  the  one-parameter  family  of  parabolas  2((+i)i  =  P  parameterized  by  p  from  0  to  1 .  The 
parabolas  have  their  vertices  colocated  at  (£,  77)  =  (  —  1,0),  i.e.,  at  the  lower  left  corner  of  the 
isosceles-like  triangle  and  their  common  tangent  the  baseline  77  =  0.  When  p=0,  the  parabola  is 
the  baseline  77  =  0.  When  p=l,  the  parabola  is  the  left  boundary  of  the  isosceles-like  triangle 
77  =  (£  a- 1)2.  Therefore,  this  one-parameter  family  of  parabolas  will  cover  the  whole  isosceles-like 
triangle  when  the  parameter  p  is  varied  from  0  to  1. 

We  need  to  find  the  shape  of  the  contours  of  the  eigenvalue  ratio  ^  =  /,  i.e., 

At  (l  +  t)  +  y/(l  +  ty-ri  , 

a2  (l  +  fwa+o2-? 

Having  set  up  the  problem  in  this  way,  it  is  easy  to  see  that  the  contours  of  the  eigenvalue  ratio 
=  /  are  the  same  one-parameter  family  of  parabolas  introduced  above  with  the  conversion 
between  the  parabola  parameter  2(f+Tp  ~  P  an^  eigenvalue  ratio  jJ-  =  /  as 


P  = 


4/ 

(*  +  l)2 


or  l  = 


2  -  p  +  i/4  -  4p 
P 


We  can  reacquaint  ourselves  with  the  previously  found  result  that  the  eigenvalue  ratio  ^  =  l 
is  1  on  the  left  boundary  of  the  isosceles-like  triangle  77  =  (f  + 1)2  and  0  on  the  baseline  77  =  0.  The 
dB-scale  contour  plot  for  this  equal  strength  arrival  case  is  presented  in  fig.  5.10  in  section  5.2. 
Because  the  contours  of  ^  =  /  pass  through  the  lower  left  corner  of  the  isosceles-like  triangle,  this 
corner  is  a  point  of  discontinuity  for  the  eigenvalue  ratio  for  the  special  case  of  equal  strength 
=  0  dB. 
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3.5  General  Case 

The  quadratic  equation  from  the  characteristic  polynomial  has  the  form 


3.5.1  Straight  Line  Contours  for  Constant  Eigenvalue  As 

From  the  quadratic  equation,  we  can  see  that  the  constant  As’  contours  remain  straight  lines  in 
the  (  and  q  plane  as  in  the  special  case  ^  =  1  as 


i+ <?)-'*’  = 


,  7Ti  .  1.  ,  .  l.TTi.  1 

2(-i)-’A  (+-(-^-)  5  + 

7T  2  Z  7T2 


The  mesh  and  contour  plots  of  Ai  are  shown  in  figures  3.7,  3.9,  and  3.11  for  —  10,20,30 
dB  respectively.  The  corresponding  plots  of  A2  are  shown  in  figures  3.8,  3.10,  and  3.12.  These 
straight  lines  have  positive  slopes  2A;(f^)~2  and  intercept  of  -(f^)-1A?  with  the  vertical  line 
passing  through  the  left  vertex 


As  increases  from  unity,  the  slope  of  the  contours  of  A2  are  smaller  as  the  range  of  A2  is  still 
between  0  and  1.  The  slope  of  the  contours  of  A!  are  larger  as  the  range  of  Ax  is  of  the  order  of 
^  rather  than  that  of  In  other  words,  the  variation  of  the  large  eigenvalue  Ai  is  essentially 
along  the  direction  of  the  horizontal  phase-dependent  variable  (,  and  the  variation  of  the  small 
eigenvalue  A2  is  essentially  along  the  direction  of  the  vertical  phase-independent  variable  7/. 


3.5.2  Parabolic  Slices  for  Constant  Phase-Dependent  Variable  (s 

For  constant  (  slices,  the  characteristic  equation  can  be  rearranged  to  that  of  a  parabola  in  A  and 
rj  as 


Figures  3.13,  3.15,  and  3.17  plot  for  A!  And  A2  vs  q  At  (=-0.50,-0.25,  0,  +0.25, +0.50  For  ^  = 
10,20,30  dB. 

We  note  that  for  given  (,  both  coordinates  of  the  vertex  of  the  parabola  increase  with  respect 
to  ,  as  their  derivatives  with  respect  to  ^  have  the  following  forms 


d 


a  [O  +  S  +  V?2^ 

1T2  IT  2 
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Figure  3.12:  Straight  line  contours  of  A2  for  ^  =  30  dB. 
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lambda  1  &  2  S  lambda  1  &  2 


Figure  3.18:  Hyperbolic  slicesrA,  and  X2  vs  £  at  t?  =  0.0,0.25,0.50,0.75,1.00  for  ^  =  30  dB. 

K2 
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Therefore,  as  increases  from  1,  those  parabolas  originally  connected  in  A  and  77  on  the 
left-half  of  the  domain  of  interest  in  the  £  and  77  plane,  i.e.,  £  <  0,  start  to  depart  from  each  other. 
Those  parabolas  on  the  right-half  of  the  domain  of  interest,  i.e.,  at  £  >  0  and  already  separated 
when  =  1,  separate  even  farther  as  ^  increases. 

The  curvature  of  each  parabola  is  2(^-)-1  at  the  vertex  which  is  outside  of  the  domain  of 
interest  as  mentioned  earlier.  For  given  ji,  the  curvature  is  independent  of  the  £  value  and  all 
the  parabolas  are  exact  translational  copies  of  each  other. 

Taking  partial  derivatives  of  the  parabolic  equation  with  respect  to  77  and  rearranging,  we  have 

a  A = _ is _ _  _  _ 

dV  2  [A  -  |(1  +  £  +  ^20]  (1+lL+fE2i)i 

1 — 4 — 77 

At  77  —  0,  we  have 

d  i 

#-A(C,0)  = - -*-?=-■ 

S1  1  +  ^  + 

For  large  we  have  the  following  approximation 

|A«'°)  -  *»• 

The  local  first-order  approximations  to  As  near  77  =  0  are  therefore 

«  A!^,0)  =  +  1  +  ,/^2£, 

772  V  1*2 

A2(^,  7?)  «  77. 

We  note  the  important  observation  that  the  smaller  eigenvalue  is  essentially  independent  of  both 
£  and  near  the  baseline  for  large  The  ratio  ^  essentially  only  affects  the  large  eigenvalue. 
The  large  eigenvalue  is  not  sensitive  to  whether  the  scenario’s  (£,77)  coordinate  is  close  to  the 
baseline,  i.e.,  whether  the  arrivals  are  highly  temporally  correlated  and/or  spatially  close  to  each 
other. 

3.5.3  Hyperbolic  Slices  for  Constant  Phase-Independent  Variable  77s 

For  constant  77  slices,  the  characteristic  equation  can  be  rearranged  to  that  of  a  hyperbola  in  A 
and  £  as 

A2“\/r2^_(1  +  7')A+~,?  =  0- 

V  7T2  X2  7T2 

In  fact,  for  each  ^  value,  this  describes  a  one  parameter  77  family  of  hyperbolas  in  A  and  £ 
that  share  common  asymptotes.  We  defer  the  derivation  and  some  of  the  characterization  of 
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these  hyperbolas  to  the  appendix.  Figures  3.14,  3.16,  and  3.18  plot  for  Ai  and  A2  vs  £  at  rj  = 
0.0,0.25,0.50,0.75,1.00  for  =  10,20,30  dB  respectively. 

For  a  given  77,  the  constant  77  slices  of  X\  and  A2  curves  as  functions  of  £  are  the  truncated  parts 
of  the  same  hyperbola  in  A  and  (.  Because  of  this  characterization,  (Ai)maj.  and  (A2 )mm  occur  at 
the  same  place,  £mar.  This  is  convenient,  as  we  are  not  yet  sure  whether  we  can  consolidate  the 
two  possible  "stress”  measures,  A2  and  into  one.  Therefore,  the  worst  point  is  always  at  the 
right  vertex  of  the  baseline,  corresponding  to  constructive  interference. 

The  higher  ^  values,  as  contrasted  to  the  special  case  of  ^  =  1,  have  wider  separation  of 
the  Ai  and  A2  surfaces.  The  truncated  parts  of  each  hyperbola  in  A  and  £  above  the  isosceles- 
right-triangle-like  region  get  closer  to  the  asymptotes,  as  the  intersection  of  the  two  asymptotes 
(-5(14- moves  to  the  left  from  (-1,0).  Above  the  baseline,  i.e.,  at  77  =  0,  the  maximum 
change  of  Ai  with  respect  to  its  mid-point  value  for  »  1  is 


4,/^ 

V  *2 

1  + 

'  IT 'i 


)• 


That  is,  the  dependence  of  As  on  £  becomes  more  insignificant  as  (p  gets  large. 


3.5.4  Contours  for  Eigenvalue  Ratio  jp 

We  need  to  know  the  shape  of  the  contours  of  the  eigenvalue  ratio  ^  =  /,  i.e., 

A,  1  -  S  ' 

where  y  assumes  the  expression  from  section  3.1 


1  - 


4(g>1 


\  [s  +  i+(g>M' 


This  is  equivalent  to 

+  1+1 
or 

4(^)7  =  41 

+  (,+1)r 

The  contours  of  the  eigenvalue  ratio  ^  =  /  for  fixed  (p  form  a  one-parameter  family  of 
parabolas  having  their  vertices  colocated  at  ((,77)  =  (-|(1  -I-  i.e.,  the  intersection 

of  the  two  asymptotes  described  in  the  last  subsection,  and  their  common  tangent  the  baseline 
77  =  0.  Compared  to  the  equipower  arrival  case,  we  note  that  this  common  vertex  has  moved  to 
the  left  from  (-1,0)  and  is  now  outside  of  the  isosceles-like  region.  As  a  parabola  has  its  largest 
curvature  at  its  apex,  the  contours  over  the  isosceles-like  region  should  have  less  curvature  as 
the  (p  increases.  This  can  be  seen  from  the  dB-scale  contour  plots  presented  in  section  5.2  in 
figures  5.12,  5.14,  and  5.16  for  power  ratio  (p  =  10,20,30  dB  respectively. 
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We  can  also  characterize  the  contours  of  the  eigenvalue  ratio  ^  =  /  for  £  =  0  from 
logic  I?  =  2l°glo(^  +  1)  -  lo§lo(^)  +  l°glo(^~^2)  *  loglo(^)  +  logl°(  (7^1)2  )• 


Therefore,  for  large  the  contours  of  the  eigenvalue  ratio  ^  =  /  for  £  =  0  as  functions  of  ^ 
and  r]  are  45°  parallel  lines  with  all  three  variables  ^  =  /,  and  tj  in  dB  scale.  This  can  be 
seen  from  the  dB-scale  contour  plot  presented  in  section  5.2  in  figure  5.32. 
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3.6  Specializing  Speiser’s  Eigenvalue  Bounds  to  Two  Sources 

As  we  have  the  exact  expressions  of  the  eigenvalues  for  .he  two-source  scenario,  we  can  examine 
quantitatively  some  general  multisource  results  by  specializing  them  to  this  case.  For  example, 
we  consider  the  singular  value  bounds  derived  by  Speiser  and  Arnold[5l]  for  the  product  of  two 
square  matrices.  Specialized  to  the  two-source  scenario  and  in  terms  of  eigenvalues  in  our  context, 
his  expression  can  be  written  as 

CT2(P)°2W  <  x2  <  min  , 

where 


T  1  4> 

\If  =  ^ 

0*  1  ’ 

and  of('l')  is  the  i-th  largest  eigenvalue  of  the  Hermitian  matrix  $ ,  and  from  section  3.1 

A,  =  <x2(P¥). 


It  is  easy  to  see  that 

This  corresponds  to  A,  for  the  special  case  £  =  0,  i.e.,  at  the  vertical  axis,  with  P  being  identity, 
the  temporally  uncorrelated  equal-powered  arrival  case.  We  also  have 


4(%)(1  ~  \P\2)  ' 

(£  +  l)2 


This  corresponds  to  A,  for  the  special  case  £  =  0,  i.e.,  at  the  vertical  axis,  with  ^  being  identity, 
the  spatially  orthogonal  arrival  case. 

Given  these  algebraic  expressions,  we  now  interpret  their  physical  significance.  In  chapter  , 
we  discussed  the  three  main  sources  for  the  signal  eigenvalue  spread  or  causes  for  small  signs 
eigenvalues.  They  are  the  high  spatial  correlation  (j>  due  to  narrow  spatial  separation  between  two 
arrivals,  the  high  temporal  correlation  p  because  of  multipaths,  and  the  high  strength  contrast 
We  note  that  0  appears  in  ty,  p  and  ^  in  P.  Speiser’s  eigenvalue  bounds  for  the  product  of  two 
matrices  tell  us  that  if  we  know  the  eigenvalues  of  $  and  P ,  then  we  can  bound  the  eigenvalues, 
A2  =  a\ (P'P),  of  the  product  matrix.  Stated  loosely,  if  we  know  the  partial  stress  from  spatial 
correlation  alone  and  that  from  temporal  correlation  together  with  the  power  contrast,  we  can 
bound  the  joint  stress  from  all  of  the  above  combined  sources.  From  the  above  expressions,  we 
know  that  the  partial  stresses  and  a\ (P)  depend  on  <t>  and  p  constrained  on  the  two  unit- 

disks  only  through  their  magnitudes,  i.e.,  radii,  but  not  the  angular  positions.  Yet  the  combined 
stresses  A2  =  o\ (P'f)  depend  also  on  the  relative  phasor  positions.  As  we  have  used  the  natural 
many-to-one  mappingfrom  ( p,<f> )  to  (Z,tj)  which  appeared  in  the  eigenvalue  formulas  in  section  3.1. 
we  note  that  the  relative  phasor  positions  appeared  only  in  the  phase-dependent  variable  £. 

For  given  \p\  and  |0|,  we  have  rj  =  (1  -  |<£|2)(1  -  ]p|2)  and  |£|max  =  \p\\4>\-  It  is  of  interest  to 
see  whether  Speiser’s  upper  and  lower  bounds  for  A2  coincide  with  A2(-|p||<?!>|, (1  -  |d>|2)(l  -  |p|2)) 
and  A2( |p||<£|, ( 1  -  |<£|2)(1  ~  W2))-  We  also  like  to  know  whether  these  bounds  are  tight  or  loose 
as  the  scenario  points  range  over  the  isosceles  region  for  the  equipowered  arrival  —  1  case  and 
the  large  power- ratio  case. 
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3.6.1  For  the  Special  Equipower  Case  =  1 
For  ^  =  1,  Speiser’s  expression  states  that 

(1  -  |p|)(l  -  H)  <  A2  <  min[(l  +  H)(l-|p|),(l-|^|)(l  +  |p|)] 

=  (i-WM)-IIpI-MI. 

These  bounds  have  exactly  the  same  expressions  for  A2  at  the  end  points  of  the  permissible  (- 
interval  between  (-|p||<£i,  (1  -  |</>|2)(1  -  |p|2))  and  (|p||<^|,(l  -  |<£|2)(1  -  |p|2))  when  =  1  as  shown 
in  the  following. 

For  a  given  77  value  and  ^  =  1,  as  (  varies  from  the  left  edge  of  the  permissible  (-interval  to 
the  midpoint  and  then  to  the  right  edge,  A2  decreases  from 

A2(-ipiM,(l-H2)(l-|p|2)) 

=  (1  -  \p\m  -  -  \p\m 2  -  a  -  mu  - 1  h2) 

=  (i  -  \m)  -  m  -  m, 

to 

A2(0,(1  -  |*|2)(1  -  |p|2))  =  1  -  yi-(i-|<AI2)(i-H2), 
to 

A2(|p||^i,(l-|^|2)(l-|p|2)) 

=  (1  +  IpIM)  -  y/(  1  +  IpIM)2  -  (1  -  mu  -  IpI2) 

=  (i  +  HH)-(H  +  H)  =  (i-M)(i-|p|). 

So  for  ^  =  1,  Speiser’s  two  bounds  for  the  equipower  case  can  be  identified  as  the  eigenvalue 
values  at  the  two  extreme  ends  of  the  permissible  phase-dependent  variable  (-interval.  This 
(interval  is  on  the  constant  phase-independent  77  horizontal  line  of  the  isosceles  region  with 
all  expressions  calculated  from  the  magnitudes  of  p  and  (j>  only.  For  given  77,  the  permissible 
(-interval  attains  its  maximum  width  when  the  (-interval  extends  from  the  left-  to  the  right- 
boundary  parabola.  The  length  of  such  maximum  (-interval  increases  with  decreasing  77,  i.e., 
when  the  baseline  is  approached.  The  bounds  for  the  small  eigenvalue  may  be  loose  when  there 
is  considerable  variation  of  the  order  of  magnitude  of  the  small  eigenvalue  along  the  maximum 
(-interval. 


3.6.2  For  the  General  Power  Ratio  2:1  Case 

*2 

For  the  general  case,  Speiser’s  expression  can  be  similarly  arranged  as 


^  A2 


Note  these  bound  expressions  involve  only  the  absolute  magnitude  of  the  two  phasors  on  the 
two  unit-disks,  i.e.,  only  the  radius  but  not  the  argument  of  the  phasors.  They  are  not  explicitly 
dependent  on  the  other  parameter 

(  =  R e(p</>*)  =  |p||<£|  cos(argp  -  arg<£), 


Z  7T2 


\4\ 


l  - 


4(^)(l-|p|2) 


( £  +  l)2 


(£  +  !)2 
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which  is  the  inner  product  of  these  two  phasors  and  consequently  contains  relative  phase  infor¬ 
mation  between  the  two.  But  as  before,  for  given  |p|  and  |<£|,  we  have  77  =  ( 1  -  |<£|2)(1  -  |p|2) 
and  |£|max  =  IpII^I-  Therefore,  as  in  the  preceding  subsection  for  the  equipower  case  ^  =  1, 
we  can  draw  the  conclusion  a  priori  that  these  bound  expressions  at  best  bound  the  lower  part 
of  the  truncated  hyperbola  in  £  and  A  over  the  permissible  phase-dependent  variable  ^-interval 
corresponding  to  the  given  value  of 


t?  =  (i-|4>|2)(i-IpI2). 

For  given  value  of  the  exact  eigenvalues  depend  on  the  two  complex  phasors  on  the  unit- 
disks,  p  and  <f>,  through  the  two  derived  bounded  real  numbers  £  and  rj.  In  contrast,  the  two 
bounds’  expressions  involve  the  magnitude  of  the  two  phasors,  i.e.,  two  different  real  numbers. 
Therefore,  it  is  much  simpler  to  browse  the  hyperbolic  slice  plots  to  visualize  the  best  non-explicitly 
^-dependent  bounds  one  can  possibly  have  than  to  analyze  the  particular  bounds  derived  by 
Speiser.  From  our  early  characterization  of  the  A  surfaces,  we  know  that  the  dependence  of  As 
on  £  becomes  more  insignificant  as  gets  large.  We  expect  that  such  ^-independent  bounds 
perform  better  in  this  region.  Naturally  we  must  bear  in  mind  that  the  original  form  of  these 
bounds  are  applicable  to  the  general  multisource  arrival  scenario  and  they  may  serve  some  useful 
purpose  there. 

We  have  used  plots  to  give  some  indication  of  the  ^-dependent  nature  of  the  eigenvalues.  We 
can  also  give  some  analytical  indication  of  the  same  phenomenon  for  the  most  important  area.  As 
the  77  =  0  baseline  is  where  A2  =  0,  we  have  more  quantitative  information  about  the  ^-dependent 
or  ^-independent  nature  of  A2  in  this  neighborhood  for  given  values  of  because  of  the  following 
analytical  characterization.  Near  77  =  0  and  at  the  three  particular  (  values,  we  have  for  ^  =  1 

•M-l,*?)  *  ooj?,  A2(0, 77)  «  A2(1,t7)%j. 

This  indicate  the  significant  ^-dependence  of  A2.  Near  77  =  0  and  for  large  we  have  instead 

A2(£,»?)  ~  V, 


i.e.,  essentially  ^-independent. 


48 


Chapter  4 


Fractional  Beamwidth  Separation 
Between  Two  Plane  Waves 
Impinging  on  a  Uniform  Line  Array 
(ULA) 


Even  though  we  only  assess  those  algorithms  applicable  to  general  sensor  array  configurations, 
we  still  need  to  understand  the  results  for  plane  waves  impinging  on  a  ULA,  especially  those 
simulation  results  that  could  be  expensive  to  duplicate  for  general  arrays.  After  all,  more  than 
half  of  the  theoretical  and  numerical  development  on  ABF  and  high- resolution  DF  are  specialized, 
occasionally  unnecessarily,  to  ULA.  Many  ULA  specific  results  can  shed  light  on  the  general  array 
as  well. 

Frequently  in  the  direction-finding  literature,  the  stress  to  a  direction-finding  scenario  in  the 
context  of  close  arrivals  is  expressed  in  terms  of  the  fractional  beamwidth  separating  the  two 
arrivals,  with  or  without  the  number  of  elements  of  the  ULA  prominently  specified.  There  seems 
to  be  the  implicit  assumption  that  the  stress  depends  mainly  on  the  fractional  beamwidth  which 
has  already  folded  into  it  the  information  about  the  number  of  sensors,  N,  used  in  the  array. 
In  other  words,  the  figures  plotted  against  fractional  beamwidth  for  a  ULA  using  some  specific 
number  of  half-wavelength-spaced  elements  are  indicative  of  the  behavior  for  another  ULA  using 
different  number  of  elements.  Therefore,  it  is  useful  to  assess  the  validity  of  this  assumption  for 
the  uncorrelated  special  case  for  ULA,  as  it  is  the  easiest  one  analyzed.  Here,  we  express  the 
stress  in  terms  of  the  smallest  eigenvalue  X2  (X  r/  =  l  -  |d>!2,  which  is  good  for  close  arrivals  and 
all  51  values  near  the  rj  =  0  baseline  region,  where  4>  is  derived  in  the  usual  fashion. 


4.1  Fractional  Beamwidth  Measurement  for  Arrival  Separation 
^frac.bw 

Let  the  elements’  location  of  the  ULA  have  the  following  coordinates, 

[(4,  Vn),  n  =  1,  •  •  • ,  N]  with  xn  =  (n  -  1  )d,  yn  =  0, 
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where  d  is  the  uniform  interelement  spacing.  Define  as  in  [23,  page  38]  or  [52,  page  48]  the  reduced 
angular  variable 

0  =  sin  a  -  sin  a0  =  2  cos  -( a  -|-  ao)  sin  -(a  -  a 0), 

where  a  and  a o  represent  the  directions  of  two  arrivals  relative  to  broadside.  Let  k  be 

2tt/ wavelength,  i.e.,  the  wave  number.  The  normalized  inner  product  of  the  two  steering  vectors 

corresponding  to  the  two  directions  of  plane  wave  arrivals  a  and  c*o  is 


<t> 


±(1  +  ei^e  +  . . .  +  e](N-\)kd6^  _ 


1  e’Nkd6  -  1 
N  eikde  -  1 


_  _j_  j±(N-i)kd6  sin (%Nkd0) 

N  sin  (^kdO) 

In  section  3.1,  we  listed  the  known  property  that  the  inner  product  of  two  steering  vectors 
for  plane  waves  impinging  onto  pairwise  symmetric  arrays  is  real.  The  above  spatial  correlation 
equation  indicates  that  the  phase  term  will  disappear  with  the  change  of  the  coordinate  origin 
from  one  end  to  the  center  of  the  ULA. 

For  a  continuous  line  array,  the  standard  beamwidth  definitions  are  inferred  from  the  first 
null  of  the  radiation  pattern  in  units  of  the  reciprocal  of  the  array  aperture.  From  the  above 
displayed  equation,  it  is  clear  that  for  the  discrete  ULA,  the  length  of  the  array  aperture  should 
be  interpreted  as  Nd,  not  ( N  —  1  )d.  For  ULA,  the  4-dB-down.bwJn-0  beamwidth  coincides  with 
the  standard  beamwidth,  see  [22,  eqns.  4.5,4.34,4.35]  for  example.  However,  the  4-dB  down 
beamwidth  definition  can  be  readily  extended  to  general  arrays  whereas  that  of  the  standard 
beamwidth  cannot.  For  example,  the  log-periodic  array  does  not  have  nulls  in  its  radiation 
pattern.  This  is  based  on  the  equivalent  Logarithmic  Period  Modulated  (LPM)  or  Hyperbolic 
Frequency  Modulated  (HFM)  waveform  result  [11]  where  ihe  sidelobe  magnitude  square  has  a 
positive  term  in  additional  to  the  familiar  sine-square  term.  We  will  also  relate  the  separation 
between  two  close  arrivals  in  terms  of  the  fractional  beamwidth  for  the  mechanical  directions  of 
arrivals  a  and  a0  to  that  in  terms  of  the  electrical  phase  angle  0. 

The  4-dB  down  beamwidth  in  the  electrical  phase  angle,  kd0,  from  element  to  element  along 
the 

ULA  (4-dB -down -bwJnJS)  is  defined  as  the  reciprocal  of  the  array  aperture  by  [5,  eqns.  6-34,11- 
12,11-14] 


4-dB -down  .bwJn  .0 


Nd 

wavelength. 


2  TT 

~Nkd ’ 


where  k  =  27r/wavelength.  The  corresponding  4-dB  down  beamwidth  relative  to  the  two  direc¬ 
tions  of  arrivals  a  and  op  ( 4-dB-down-bw-rel.to-(a  +  ato)/2)  is  defined  as  the  reciprocal  of  the 
effective  array  aperture  by 


4-dB -down -bw-reLto.(a  +  ao)/2 


Ndcos  l(a  +  ao) 

wavelength 

2tt 

Nkd cos  i(a  +  ao) 


The  separation  between  two  arrivals  can  be  expressed  as  0]racJ>w  in  terms  of  fractional 
beamwidth  by 


0 


frac.bw 


_ 0 _ 

4-dB.down.bwJn-0 


NkdO  _  Nd0 
2x  wavelength’ 
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or  equivalently 


2  T 

—  ~N  Qfrac.bw 


The  terminology  makes  sense  as  for  close  arrivals 


a  -  o0 


4-dB-down-bw-reLto _  (a  +  ao)/2 


2  sin  |(a  -  ao) 

Nd  cog  J-(a+oro)1  1 

•avelength 

2 sin  ^(a  -  a0) cos  ^(a  +  a0) 


Nd 


L wavelength] 


=  (sin a- sin qq) 


Nd 


[wavelength] 

=  0  x  apertureJn.units.of-wavelength  =  0frac  ]~)W- 

The  above  development  is  similar  to  that  reported  by  Evans,  Johnson,  and  Sun  [15]. 

The  normalized  inner  product  of  the  two  steering  vectors  corresponding  to  the  two  directions 
of  arrivals  a  and  ao  is 

<f>  =  —  ejn^efrac.bw  J^6frac.bw) 

N  s™(-h*0frac_bw) 

Therefore,  we  have 

s\ti(k6  fracbw) 
s,n(Nw^frac.bw)  I 

For  closely  spaced  arrivals,  the  above  equation  can  be  approximated  by 


\<t>\  =  — 
m  N 


\<}>\  « 


*  6  i^^frac-bw]  + 


1 


i  r  i 

~  6  N^^frac-bw  + 


2  K  jy^X^fracJjw)  • 


For  very  close  arrivals,  <f>  «  1,  i.e.,  i'  location  on  the  unit-disk  is  on  the  circumference.  The 
horizontal  axis  has  the  following  intei,  station  now:  £  «  Re(p)  =  |p|  cos(argp). 


4.2  The  Smaller  Eigenvalue  A2’s  Negligible  Dependence  on  the 
Number  of  Elements  N  after  Using  0frac  for  the 
Equipower  Uncorrelated  Case 


For  uncorrelated  and  equal  strength  arrivals,  the  only  stress  to  the  direction-finding  scenario  is 
the  angular  separation  between  the  arrivals,  therefore 


a2  =  i  -  \<f>\  =  l  -  I 


sin(7rfy-rac_fow) 

s'^N^frac.bw) 


For  this  case,  a  discussion  of  the  smaller  eigenvalue  A2  directly  results  from  that  for  the  magnitude 
\4>\  of  the  generalized  cosine  of  the  steering  vectors  of  the  two  arrivals.  An  approximation  for  A2 
for  closely  spaced  arrivals  can  also  be  obtained  via 


A2 


1  ~  |fla 
2 


1  sin(xflfrac_fow) 

Asin(^*0frac_bw) 


jy2 X^fyjrac-bw)  • 
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All  these  expressions  for  <f>  and  t)  depend  on  both  9frac  5W  and  N. 

For  small  0fl.ac  jjW  associated  with  the  stressing  scenario  of  close  arrivals,  say  0frac  5W  <  0.1, 
the  last  equation  shows  that  A2  increases  quadratically  with  #frac  5W  and  increases  with  increasing 
N.  See  figure  4.1.  The  dependence  on  N  is  noticeable  only  for  small  N  but  becomes  rapidly 
negligible  for  large  N.  See  figure  4.2.  The  last  equation  indicates  that  the  variation  of  &frac 
for  N  >  10  is  only  1%.  The  A2  value  at  N  =  3  is  90%  of  that  at  IV  =  00.  Therefore,  one  can 
safely  ignore  the  A2’s  additional  dependence  on  N. 

As  6frac  ijW  increases,  say  0.1  <  #frac  <  1.0,  the  unbounded  parabola  approximation  is 
obviously  less  valid  as  the  (1  -  sinc2)-like  function  is  bounded  above  by  1.  We  note  that  the 
stress  in  terms  of  A2  is  among  the  lowest  possible  near  &{rac  =  1.0.  When  the  two  arrivals  are 
separated  by  exact  multiples  of  the  Raleigh  beamwidth,  the  inner  product  between  the  two  steering 
vectors  is  zero,  i.e.,  the  two  steering  vectors  are  orthogonal.  The  A2’s  additional  dependence  on 
N  after  using  9frac  is  even  smaller  than  the  previous  6frac  bw  —  0.1  case.  See  figures  4.3  and 
4.4.  We  note  that  in  both  figures,  A2  =  1  at  0frac bw  =  1.0  regardless  of  the  values  of  N,  as  it 
should. 

If  we  continue  the  plot  beyond  0frac  =  1.0,  the  ripple  phenomena  of  the  A2  versus  #frac  j,w 
will  show  up  prominently  as  indicated  in  figures  4.5  and  4.6.  This  is  associated  with  the  sidelobes 
of  ULA,  which  depends  on  both  9frac  bw  and  N.  That  is,  the  (£,  rj)  will  stay  near  (0,1),  the 
apex  of  the  isosceles-right-triangle-like  region.  The  stress  to  such  a  direction-finding  scenario 
in  terms  of  A2  can  only  come  from  the  two  arrivals’  strength  ratio  which  will  show  up 

in  the  A1/A2  ratio  and/or  source  amplitude  correlation  p  such  as  from  multipath.  Overlooking 
the  ripples,  the  (1  -  sinc2)-like  function  gets  closer  to  1  as  0frac  bw  increases.  Therefore,  the 
dependence  on  N  is  not  significant  either.  We  also  notice  the  grating  lobe  phenomena  such  as  at 
(9frac  =  (3, 3), (6, 3), (9, 3), (6, 6),  (9,9)  in  figures  4.7  and  4.8.  See  [52]  for  explanations. 

In  figure  4.8,  we  have  chosen  to  plot  the  A2  =  1  -  |<£|  slices  vs  N  at  0frac  =  0.5  :  1  :  9.5,  which 
correspond  to  the  sidelobes  of  the  CBF. 

However,  because  of  the  ripple  effect,  the  smaller  eigenvalue  may  even  decrease  when  more 
sensors  are  used  while  everything  else  is  the  same.  This  will  cause  some  undesirable  phenom¬ 
ena  when  noise  is  brought  in  and  a  unit-weighting  signal  space  method  is  used.  For  moderate 
SNRs,  the  error  variance  of  each  eigenvector  calculated  from  the  sample  data  covariance  matrix  is 
roughly  proportional  to  that  particular  eigenvalue.  Ottersten,  in  fact,  produced  some  interesting 
asymptotic  results  illustrating  the  detrimental  effect  to  the  DOA  parameter  standard  deviation 
when  plotted  as  a  function  of  the  number  of  sensor  elements  for  the  unity-weighting  algorithm  for 
ULA[38,  pages  68-69].  As  he  fixed  the  mechanical  angle  separation  between  the  two  arrivals,  the 
fractional  beamwidth  changes  as  the  number  of  sensor  elements  varies.  The  (0frac  bWiN)  points 
trace  out  a  curve  in  the  base  plane  of  the  mesh  and  contour  plots  in  figures  4.5  and  4.6. 

In  conclusion,  the  fractional  beamwidth  notion  is  a  useful  concept  for  assessing  the  stress  to 
a  direction-finding  scenario  in  terms  of  the  normalized  smaller  eigenvalue  A2.  It  consolidates  the 
effects  of  arrival  separation,  sensor  element  spacing,  wavelength,  and  number  of  sensor  elements 
into  one  entity. 
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Figure  4.6:  Contour  plot  of  X2  =  1  —  \<t>\  over  1  <  0frac  bw  <  10  and  N=3:l:22  for  ULA. 
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Chapter  5 

Plots  of  Signal  Eigenvalues  and  their 
Ratio  in  dB  Scale 


The  plots  of  the  same  information,  or  function,  using  direct  scale  and  dB  scale  in  the  vertical 
axis  have  their  individual  merits.  The  direct-scale  versions  show  the  parabolic  and  hyperbolic 
sections  as  well  as  the  straight  line  contours  in  their  natural  coordinates  as  derived  from  the 
analytical  studies.  Versions  using  the  dB  scale  better  illustrate  the  multiplicative  dependence 
of  the  eigenva’ues  on  some  of  the  three  independent  parameters.  In  addition,  the  effect  of  the 
smallest  dominant  eigenvalue  is  relative  to  the  threshold  setting,  which  is  frequently  expressed  in 
dB  scale;  the  condition  number  can  be  exhibited  more  compactly  in  this  way. 

While  the  parabola  sections  are  all  translation  copies  of  each  other  for  given  ^  values  in 
the  direct-scale  coordinates,  one  cannot  recognize  such  relations  from  the  dB-scale  plots.  For 
example,  the  vertical  intervals  (0,9)  and  (1,10)  are  translation  copies  of  each  other  and  can  be 
recognized  as  such  easily  in  the  direct-scale  plot.  But  their  images  in  the  dB  scale  are  qualitatively 
different  because  of  the  nonlinear  nature  of  the  logarithmic  tunction,  especially  near  zero,  which 
is  so  important  for  our  study. 

In  section  5.1,  we  first  display  the  slices  of  the  two  signal  eigenvalues  Aj  and  A2  in  dB  along 
constant  phase-independent  variable  77  and  along  constant  phase-dependent  variable  £  for  power 
ratio  ^  =  0,10,20,30  dB.  Next,  their  ratio,  which  is  the  difference  in  dB  of  the  large  and  small 
signal  eigenvalues,  will  be  examined.  In  sections  5.2  and  5.3,  we  will  find  that  the  variation 
of  signal  eigenvalue  ratio  ^  across  the  phase-dependent  variable  £  is  not  significant.  We  focus 
attention  to  the  condition  number’s  behavior  over  the  vertical  axis  £  =  0  of  the  isosceles  region 
in  section  5.4. 

While  the  parabola  and  hyperbola  characterization  described  in  sections  3.4  and  3.5  pertains 
only  to  linear-scale  plots,  we  retain  for  convenience  such  terminology  even  for  those  functional 
representation  in  the  semilog  coordinates  used  in  this  chapter. 

As  77  approaches  0,  the  small  eigenvalue  A2  and  the  eigenvalue  ratio  ^  tend  toward  -  and 
+  00  dB  respectively.  Similarly,  as  77  approaches  0  for  the  equipower  case,  ^  =  0  dB,  both  the 
upper  Aj  and  lower  A2  hyperbolas  tend  toward  their  asymptotes  intersecting  at  -00  dB.  Because 
we  cannot  display  ±00  dB,  we  choose  to  stop  at  a  small  value  of  7/  =  -13  dB,  which  corresponds 
to  spatially  orthogonal  arrivals  with  temporal  correlation  of  97.5%,  or  temporally  uncorrelated 
arrivals  with  a  |  fractional  beamwidth  spacing  for  a  ULA. 
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5.1  Slices  of  Signal  Eigenvalues  Ai  and  A2  in  dB  along  Constant 
Phase-Independent  Variable  rj  and  along  Constant  Phase- 
Dependent  Variable  £  for  Power  Ratio  ^  =  0, 10, 20, 30  dB 

We  plot  the  two  signal  eigenvalues  Ai  and  A2  as  functions  of  the  phase-dependent  variable  f  and 
the  phase-independent  variable  77  for  four  given  power  ratio  ^-s. 

•  Ai  And  A2  At  £  =  -0.5, —0.25,0, +0.25, +0.5  For  Four  ^  in  figures  5.1,  5.2,  5.3,  and  5.4. 

•  Ai  And  A2  At  77  —  0.05,0.25,0.50,0.75, 1.00  For  Four  ^  in  figures  5.5,  5.6,  5.7,  and  5.8. 

For  ^  =  0  dB,  the  upper  and  lower  sheets  of  hyperbolas  and  parabolas,  the  latter  with 
negative  £,  are  all  contiguous.  They  do  not  appear  so  only  because  of  the  coarse  grid  systems  we 
used  in  making  these  plots,  i.e.,  the  left  boundary  curve  of  the  isosceles  region  does  not  fall  on 
the  grid  points.  We  mentioned  at  the  beginning  of  this  chapter  that  it  is  not  possible  to  display 
the  asymptotes  simultaneously  because  they  intersect  at  -00  dB.  We  note  that  straight  lines  in 
the  linear-scale  plot  become  curved  in  the  semilog  plot.  Over  the  apex  of  the  isosceles  triangle 
both  Aj  and  A2  show  the  common  0-dB  value  as  they  should,  in  both  hyperbolic  and  parabolic 
slice  plots. 

When  increases  from  0  dB,  over  the  apex  of  the  isosceles  triangle  Ai  and  A2  are  separated 
farther  and  farther  duplicating  this  increase.  Because  of  the  normalization  convention  adopted  for 
the  eigenvalues  with  respect  to  the  sensor  level  signal  power  of  the  weaker  source  in  section  3.1, 
over  the  apex  of  the  isosceles  triangle  A2  is  always  at  0  dB.  The  large  eigenvalue  A!  curves  are 
bundled  closer  in  both  directions  of  the  phase-independent  variable  77  and  the  phase-dependent 
variable  £.  The  variation  of  the  small  eigenvalue  A2  appears  more  and  more  only  in  the  phase- 
independent  variable  77.  As  the  eigenvalue  ratio  in  dB  is  the  difference  between  the  large  and 
the  small  eigenvalues  in  dB,  the  variation  in  the  ratio  will  be  essentially  from  that  of  the  small 
eigenvalue  in  the  phase-independent  variable  77,  as  will  be  seen  from  the  next  three  sections. 
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Figure  5.4:  Ai  and  A2  in  dB  vs  r/  at  £  =  -0.5, -0.25,0,  +0.25,  +0.5  for  ^  =  30  dB. 
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Figure  5.5:  Ai  and  X2  in  dB  vs  £  at  r}  =  0.05,0.25,0.50,0.75,1.00  for  -*■  =  0  dB. 
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Figure  5.6:  A|  and  A2  in  dB  vs  £  at  ?/  =  0.05,0.25,0.50,0.75,  1.00  for  ^  =  10  dB. 
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Figure  5.7:  Aj  and  A2  in  dB  vs  f  at  rj  =  0.05,0.25,0.50,0.75, 1.00  for  ~i  =  20  dB. 
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Figure  5.K:  A,  and  A2  in  dB  vs  £  at  r;  =  0.05, 0.25, 0.50. 0.75,  1 .00  for  =  30  dB. 
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5.2  Four  Sets  of  Signal  Eigenvalue  Ratio  Plots  for  Power  Ratio 

^  =  0,10,20,30  dB 

7T 2  7  7  7 

We  provide,  for  each  of  the  four  given  power  ratio  ^ ,  mesh  and  contour  plots  of  signal  eigenvalue 
ratio  as  function  of  the  phase-dependent  variable  £  and  the  phase-independent  variable  r 7 
over  the  rectangular  grid  [-l:l:0.05]x[0.05:l:0.05],  next  slice  plots  along  constant  phase-dependent 
variable  £  and  then  along  constant  phase-independent  variable  77.  The  vector  [-1:1:0.05]  contains 
the  numbers  from  -1  to  1  with  0.05  as  increment  in  £.  The  vector  [0.05:1:0.05]  contains  the 
numbers  from  0.05  to  1  with  0.05  as  increment  in  77. 

•  Mesh  and  contour  plots  of  in  dB  as  function  of  £  and  77  in  figures  5.9,  5.11,  5.13,  and 
5.15,  and  in  figures  5.10,  5.12,  5.14,  and  5.16  respectively. 

•  Slices  of  ^  in  dB  at  f  =  —0.5, -0.25,0, +0.25, +0.5  for  ^  =  dB  in  figures  5.17,  5.19,  5.21, 
and  5.23. 

•  Slices  of  ^  in  dB  at  ?7=0.05,0.25,0.50,0.75,1.00  for  =  dB  in  figures  5.18,  5.20,  5.22,  and 
5.24. 

For  =  0  dB,  the  mesh  plot  of  the  signal  eigenvalue  ratio  in  figure  5.9  indicates  that 
the  ratio  over  the  left  boundary  curve  of  the  isosceles  region  is  0  dB,  as  it  should  be.  We  expect 
the  slices  of  the  eigenvalue  ratio  ^  surface  for  negative  constant  £,  phase- dependent  variable,  in 
figure  5.17,  and  for  all  constant  77,  phase-independent  variable,  in  figure  5.18,  to  intersect  the  0-dB 
abscissa  axis.  Again,  they  do  not  appear  to  be  so  only  because  of  the  coarse  grid  systems  used  in 
making  these  plots,  i.e.,  the  left  boundary  curve  of  the  isosceles  region  does  not  fall  on  the  grid 
points. 

We  notice  from  these  two  slice  plots,  figures  5.17  and  5.18,  that  the  eigenvalue  ratio  ^  tends 
toward  00  dB  as  77  approaches  0.  This  illustrates  the  00  discontinuity  of  the  eigenvalue  ratio  ^ 
at  the  lower  left  corner  of  the  isosceles  region.  The  limit  approaching  from  the  baseline  of  the 
isosceles  region  is  00,  while  that  approaching  from  the  left  boundary  curve  of  the  isosceles  region 
is  0.  Note  at  this  corner,  both  signal  eigenvalues  Ai  and  A2  axe  zero  (-00  dB)  and  their  ratio 
is  indeterminate.  This  may  be  related  to  the  difficulty  mentioned  in  subsection  2.3.3.  We  have 
discussed  in  subsection  3.4.4  that  the  contours  of  the  eigenvalue  ratio  ^  =  /,  in  figure  5.12,  are 
the  same  one-parameter  family  of  parabolas  introduced  there  with  the  conversion  between  the 
parabola  parameter  2'({+i)2  =  P  an^  eigenvalue  ratio  =  l  as 


Because  of  the  difficulty  of  displaying  00  and  discontinuity,  the  apparent  absence  of  these  at  the 
left  corner  in  figures  5.9,  5.17,  and  5.18  should  be  interpreted  with  care. 

When  gets  larger,  as  in  figures  5.10,  5.12,  5.14,  and  5.16,  the  contours  get  more  and 
more  horizontal,  which  means  that  the  dependence  of  the  signal  eigenvalue  ratio  ^  on  the  two 
parameters  £  and  77  will  be  mainly  from  the  phase-independent  variable  77.  In  subsection  3.4.4,  we 
found  that  the  contours  of  the  eigenvalue  ratio  ^  =  l  form  a  one-parameter  family  of  parabolas 

having  their  vertices  colocated  at  (£,77)  =  (-|(1  +  ^)/ •>«*  figures  5.10,  5.12,  5.14,  and 
5.16.  Furthermore,  in  the  3-d  (£,77,^)  coordinates,  these  colocated  vertices  of  the  parabola 
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contour  slices  when  projected  onto  the  coordinate  plane  coincide  on  the  shared  £-axis  with 
the  intersection  of  the  common  asymptotes  for  the  77  parameter  family  of  hyperbola  slices  in  X 
and  £  for  constant  phase- independent  variable  77s  when  projected  onto  the  (£,A)  coordinate  plane. 
These  families  of  hyperbola  slices  are  displayed  in  figures  3.6,  3.14,  3.16,  and  3.18. 

The  comments  at  the  end  of  section  5.1  that  the  variation  in  the  ratio  is  essentially  from  that 
of  the  small  eigenvalue  in  the  phase-independent  variable  77  can  be  confirmed  by  comparing  the 
ratio  plots  in  figures  5.18,  5.20,  5.22,  and  5.24  with  the  small  eigenvalue  plots  figures  5.5,  5.6,  5.7, 
and  5.8. 

From  figures  5.9,  5.11,  5.13,  and  5.15,  one  observes  that  the  mesh  surfaces  are  lifted  away  from 
the  0-dB  plane  and  appear  flatter  and  flatter  when  ^  gets  larger.  This  must  be  interpreted  with 
care  because  of  the  inadequate  display  coverage  with  77  stopped  at  -13  dB  and  the  perspective 
angle  showing  the  surface  at  ^  =  30  dB  flatter  than  it  is.  While  there  is  much  less  curvature  as 
compared  to  the  equipower  arrival  ^  =  0  dB  case,  there  remains  considerable  variation  of  the 
eigenvalue  ratio  along  constant  phase-independent  variable  77.  These  plots  are  produced  by  the 
PC-MATLAB  software  which  auto-scales  the  data  to  fit  a  plotting  area  of  given  size.  If  a  common 
absolute  scale  is  used  for  all  these  mesh  plots  for  power  ratio  ^  =  0,10,20,30  dB,  the  last  three 
mesh  surfaces  would  appear  more  as  translational  copies  of  each  other.  This  can  be  seen  from 
examining  its  accompanying  slice  plots  along  constant  phase-dependent  variable  £  and  then  along 
constant  phase-independent  variable  77.  There,  the  vertical  scales  are  more  comparable  as  the  0 
dB  reference  is  not  included  as  in  the  mesh  plots.  These  will  be  made  more  clear  in  section  5.4. 

For  all  slice  plots,  we  note  that  the  apex  of  the  isosceles  region  always  occupies  the  lower  right 
hand  corner  for  the  upper  slice  plots  of  constant  phase-dependent  variable  f,  as  in  figures  5.17, 
5.19,  5.21,  and  5.23,  and  the  middle  of  the  abscissa  axis  for  the  lower  slice  plots  of  constant  phase- 
independent  variable  77  as  in  figures  5.18,  5.20,  5.22,  and  5.24.  At  the  apex,  the  signal  eigenvalue 
ratio  always  equal  to  the  power  ratio  as  it  should. 
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Figure  5.13:  Mesh  plot  of  ^  in  dB 


Figure  5.14:  Parabolic  in  dB  as  contour  pi 
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Figure  5.17:  Slices  of  ^  in  dB  at  £  =  -0.5, -0.25,0, +0.25, +0.5  for  ^  =  0  dB. 


Figure  5.18:  Slices  of  ^  in  dB  at  t?=0.05, 0.25,0.50, 0.75,1.00  for  ^  =  0  dB. 
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Figure  5.19:  Slices  of  ^  in  dB  at  £  =  -0.5, -0.25, 0, +0.25, +0.5  for  a  =  10  <1B. 


Figure  5.20:  Slices  of  ^  in  dll  at  t/=0. 05,0.25.0. 50, 0.75.1 .00  for  ^  =  10  dB. 
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Figure  5.21:  Slices  of  Q  in  dB  at  £  =  -0.5,-0.25,0,4-0.25,4-0.5  for  £l  -  20  dB. 


Figure  5.22:  Slices  of  ^  in  dB  at  77=0.05,0.25,0.50,0.75,1.00  for  =  20  dB. 
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5.3  Plots  of  Signal  Eigenvalue  Ratio  ^  Over  a  Small  Phase- 
Independent  Variable  rj  =  0.05  as  Function  of 
Phase-Dependent  Variable  £  and  Power  Ratio  ^(dB) 

We  plot,  for  a  small  constant  77  =  0.05,  signal  eigenvalue  ratio  as  function  of  the  phase- 
dependent  variable  £  and  power  ratio  then  plot  slices  along  constant  power  ratio  ^  and  then 
along  constant  phase- dependent  variable  £. 

•  ^  at  ^(dB)=0,5,10,15,20,25,30  for  77=0.05  in  figure  5.25. 

•  at  £  =  —0.50,  —0.25,0,  +0.25,  + 0.50  for  77=0. 05  in  figure  5.26. 

•  Mesh  and  contour  plots  of  ^  in  dB  for  77=0.05  in  figures  5.27,  and  5.28. 

As  we  mentioned  at  the  beginning  of  the  chapter,  we  choose  to  stop  at  77  =  -13  dB.  From 
the  slice  plot  of  the  signal  eigenvalue  ratio  ^  along  constant  power  ratio  jn  figure  5.25,  we  see 
that  the  spread  in  the  eigenvalue  ratio  is  smaller  at  the  constructive  interference  side,  £  >  0,  fhan 
at  the  destructive  interference  side,  £  <  0.  The  eigenvalue  ratio  is  an  increasing  function  of  the 
phase-dependent  variable  f  but  this  dependence  levels  off  as  the  power  ratio  ^  increases. 

From  the  slice  plot  of  the  signal  eigenvalue  ratio  ^  along  constant  phase-dependent  variable 
£  in  figure  5.26,  we  see  that  the  spread  in  eigenvalue  ratio  narrows  quickly  as  the  power  ratio 
increases.  This  means  that  the  dependence  on  the  phase-dependent  variable  £  becomes  negligible. 
The  trend  is  a  45°  line  in  the  log-log  plot,  which  means  that  the  signal  eigenvalue  ratio  ^  becomes 
directly  proportional  to  the  power  ratio  ^  in  the  direct  scale. 

The  above  observations  are  reinforced  in  the  contour  plot  of  figure  5.28  showing  signal  eigen¬ 
value  ratio  ^  as  a  function  of  the  phase-dependent  variable  f  and  power  ratio  The  variation  of 
the  contour  is  mainly  in  the  power  ratio  ^  direction.  Its  change  in  the  phase-dependent  variable 
£  direction  is  pronounced  only  for  ^  <  15  dB. 

Except  when  a  paper  focuses  on  the  effect  of  the  temporal  phase  difference,  the  phase  of  the 
source  correlation  is  frequently  implicitly  assumed  to  be  zero.  Often  we  have  to  infer  this  from  the 
following  clues:  the  source  correlation  coefficient  used  in  evaluating  asymptotic  performance  is  a 
positive  real  (without  discussion),  or  identical  complex  source  amplitudes  for  the  two  arrivals  are 
used  in  numerical  simulation.  The  position  of  a  stressful  scenario  point  in  the  isosceles  region  can 
be  anywhere  on  the  right-half  of  the  constant  77  line  parallel  to  the  close-by  baseline.  For  these 
highly  correlated  equipower  arrivals,  we  will  show  that  the  effect  on  the  eigenvalue  ratio  from 
changing  the  phase  of  the  source  correlation  from  90°  to  0°  can  be  as  high  as  a  6-dB  increase. 
The  following  paragraphs  explain  the  possible  6-dB  increase  of  the  eigenvalue  ratio  adjustment. 

The  large  eigenvalue  Aj  near  the  baseline  corresponding  to  temporally  highly  correlated  arrivals 
can  be  estimated  using  the  asymptotes’  behavior,  because  the  cross  sections  of  the  As’  surfaces  for 
constant  77s  are  hyperbolic  slices,  as  derived  in  subsection  3.1.3  and  shown  in  figure  3.6.  Depending 
on  the  relative  location  of  the  scenario  point  to  the  right  boundary  curve  of  the  isosceles  region, 
the  effect  on  the  large  eipmvalue  Aj  from  changing  the  phase  of  the  source  correlation  from  90° 
to  0°  can  be  as  high  as  a  3-dB  increase. 

The  small  eigenvalue  A2  can  be  estimated  from  the  following  fact  derived  in  subsection  3.4.2 
and  shown  in  figure  3.5.  Near  77  =  0  and  at  the  three  particular  £  values,  we  have 

A2(-M)  ~  0077,  A2(0,77)ss^,  A2(l, 77)  « 


73 


HP  UI  z  BpquiBi/1  cpqaiBi 


Depending  on  the  relative  location  of  the  scenario  point  to  the  right  boundary  curve  of  the  isosceles 
region,  the  effect  on  the  small  eigenvalue  A2  from  changing  the  phase  of  the  source  correlation 
from  90°  to  0°  can  be  as  high  as  a  3-dB  decrease.  See  figure  3.6  to  reinforce  the  above  discussion. 

For  the  case  of  almost  orthogonal  spatial  arrivals,  the  scenario  point  is  near  the  vertical  axis 
of  the  isosceles  region  to  start  with,  no  adjustment  of  the  eigenvalue  ratio  is  needed  as  we  move  to 
the  0°  temporal  correlation.  When  a  scenario  is  stressed  both  temporally  and  spatially,  then  the 
relative  position  of  the  (£,  77)  on  the  constant  77  line  as  measured  by  the  fraction  distance  between 
the  vertical  axis  to  the  right  boundary  of  the  isosceles  region  is  almost  one.  In  this  case,  the  full 
6-dB  adjustment  is  needed. 

This  can  be  seen  from  the  following.  Given  the  complex  pair  (4>,  p),  we  calculate  the  real  pair 
(£,77)  from 

7?  =  (1  -  |0i2)(l  -  H2), 

£  =  R e(p<t>m)  =  |p||<£j  cos(argp  -  arg<£). 

For  this  77,  the  boundary  point  is  defined  in  section  3.2  by 

|£|max  =  (jHl0|)max  =  1  —  \/V- 

The  fractional  distance  for  a  0°  temporal  correlation  scenario  point  is  expressed  as 

Hi  \m 

l£lmax  given  v  ^  ~  \/(l  -  |0|2)(1  ~  IpI2) 

When  both  \p\  ss  1  and  \<fr\  «  1,  we  have 

— - — - ~  \p\\4>\(i  +  2^/(1  -  |<£|)(l  -  |p|) )  ~  \p\\<f>\  %  l. 

1^1  max  given  r, 

The  above  discussions  about  the  effect  of  6-dB  increase  on  the  eigenvalue  ratio  for  highly- 
correlated  equipower  arrivals  from  changing  the  phase  of  the  source  correlation  from  90°  to  0° 
can  be  correlated  with  figures  5.25  and  5.26  for  the  small  phase-independent  variable  ij  =  0.05. 


76 


5.4  Plots  of  Signal  Eigenvalue  Ratio  ^  Over  the  Vertical  Axis 
£  =  0  as  Function  of  Phase-Independent  Variable  r]  and 
Power  Ratio  £(dB) 

Having  ascertained  that  for  ^  >  10  dB  the  variation  of  signal  eigenvalue  ratio  y^  across  the 
phase-dependent  variable  £  is  not  significant,  we  can  focus  attention  to  the  condition  number’s 
behavior  over  the  vertical  axis  f  =  0  of  the  isosceles  region. 

We  plot  for  £  =  0,  signal  eigenvalue  ratio  ^  as  function  of  the  phase-independent  variable 
Tj  and  power  ratio  next  plot  slices  along  constant  power  ratio  ^  and  then  along  constant 
phase-independent  variable  tj. 

There  are  two  versions  of  these  plots.  The  first  version  has  the  phace-independent  variable  rj  in 
dB  scale  covering  [-40:0:1],  following  the  convention  explained  in  the  previous  subsection.  The  sec¬ 
ond  version  has  the  phase-independent  variable  rj  in  direct  or  linear  scale  covering  [0.05:1.00:0.05]. 


5.4.1  With  Phase-Independent  Variable  tj  in  dB  Scale 

Signal  Eigenvalue  Ratio  y^(dB)  vs  Phase-Independent  Variable  77(dB)  Curves  at  7 
Constant  Power  Ratio  £(dB)for£  =  0 

From  figure  5.29  for  £  =  0  with  seven  curves  of  y'-(dB)  vs  7j(dB)=-40:l:0 
at  y£(dB)=0,5,10,15,20,25,30,  we  ha\e  the  following  observations. 

At  the  0-dB  right  edge  with  t?=1  which  corresponds  to  the  apex  of  the  isosceles-right-triangle- 
like  region,  each  of  the  seven  curves  corresponding  to  the  7  ^(dB)  sections  intersects  the  y^  dB 
axis  at  exactly  the  same  values  of  each  ^  curve’s  labelling,  as  they  should.  From  tj  <  -10  dB  or  < 
0.10  down,  i.e.,  except  for  the  last  two  vertical  grid  intervals,  these  7  curves  are  almost  all  straight 
lines  with  slope  of  -1  in  this  log-log  plot. 

We  notice  the  top  3  curves  corresponding  to  y^=30,  25,  and  20  dB  passing  through  the 
intersections  of  horizontal  and  vertical  grid  lines.  For  example,  the  top  curve  with  y^=30  dB 
passes  through  these  grid  points:  (tj,  -  (0,30),  (-10, 40),  (-20, 50),  (-30, 60),  and  (-40,70). 
Note  the  equal  spacing  among  the  three  curves.  When  ^  >  10  dB,  we  have 


^2  7f2 


or  that  the  stress  from  tj  and  ^  is  multiplicative  (divisive)  in  the  direct  scale.  In  dB  scale,  this 
is  additive,  i.e., 


—-(in  dB)  ss  (  — )(in  dB)  -  77(in  dB). 

A2  7T  2 

This  approximation  starts  to  degrade  when  ^  is  small  and  rj  is  close  to  1.  So  if  this  additive 
or  multiplicative  relation  were  perfect  instead  of  being  just  approximate,  the  ^-  =  10-dB  curve  will 
intersect  at  the  left  vertical  line,  tj  =  -40  dB,  to  arrive  at  y^=50  dB,  instead  of  the  51-dB  point. 
But  we  can  see  that  the  approximation  even  for  the  j^  =  10-dB  curve  is  excellent. 

For  tj  larger  than  0.10,  the  lowest  2  curves  corresponding  to  y^=0  dB  and  5  dB  exhibit 
noticeable  curvature,  especially  for  f^=0  dB  near  tj~  1.  An  original  10-dB  separation  at  tj  =  0  dB 
between  the  two  curves  y1  =  0  dB  and  y1  =  10  dB  is  reduced  to  5  dB  for  tj  <  -5  dB. 


lambda  1 /lambda  2  in  dB  lambda  1/Iambda_2  in  dB 


Figure  5.29:  ^(dB)  vs  r/(dB)=-40:l:0  at  £(dB)=0,5,10,15,20,25,30  for  £=0. 


pi_l/pi_2  in  dB  at  eta  =  cnnstant  sections 


Figaro  5.30:  ^(dB)  vs  ^(dB)  =  0:l  :30  at  r/(dB)-  —10. -30. -20. -10.0  for  (=0. 

7R 


For  the  ^-=0-dB  curve  near  77=1,  the  apex  point  area,  the  stress  to  the  direction-finding 

scenario  due  to  multipath  correlation  between  arrivals  and/or  closeness  in  spatial  separation 

increases  much  faster  than  large  dB  curves  as  77  decreases.  This  is  caused  by  the  infinite  slope 

of  A  with  respect  to  77  at  the  apex  for  the  equal-power  arrival  case.  Because  of  this  and  the 

variation  of  in  the  £  direction,  we  need  to  examine  the  surfaces  of  Ai  or  A2  or  their  ratio  in 

detail  for  the  rL=0-dB  case. 

"  2 

When  21.  «  0  dB  and  77  is  less  than  0.10  (or  -10  dB),  we  observe  that 

-i(in  dB)  as  — 77(in  dB)  -f  6. 

A2 

This  extra  6  dB  accounts  for  the  presence  of  the  infinite  slope  of  A  with  respect  to  77  at  the 
apex  for  the  equal-power  arrival  case. 

Recall  that  these  were  predicted  in  subsection  3.3.3. 

Signal  Eigenvalue  Ratio  j^(dB)  vs  Power  Ratio  ^-(dB)  Curves  at  Five  Constant 
Phase-Independent  Variable  ?7(dB)  for  £  —  0 

From  figure  5.30  for  £  =  0  with  five  curves  of  ^-(dB)  vs  ^(dB)=0:l:30 

at  77(dB)=  —40,-30,-20,-10,0,  some  of  the  just  mentioned  observations  can  be  made  clearer 
from  a  different  perspective. 

The  lowest  curve  corresponds  to  77  =  0  dB,  i.e.,  at  the  apex  point.  This  perfect  straight  line 
with  slope  of  1  indicates  the  identical  relation  between  ^  and  as  it  should.  The  five  constant 
77  curves  are  equally  spaced  in  dB.  The  intersections  of  these  curves  with  the  vertical  grid  lines 
reflect  how  well  the  corresponding  equal  spacing  can  be  carried  over  for  ^  for  different  77  values. 
We  see  that  equal  spacing  is  maintained  much  better  for  ^  >  10  dB  than  for  <  10  dB.  In 
particular,  near  the  left  vertical  edge,  which  corresponds  to  equal  strength  arrivals,  we  can  see 
that  the  stress  rises  much  faster  when  77  decreases  from  0  dB,  as  the  DOAs  get  closer  and/or 
source  amplitudes  get  more  correlated.  On  the  other  end,  the  rate  of  stress  rising  with  respect  to 
the  ^  abscissa  near  this  left  boundary  is  milder  than  other  ^  values.  In  other  words,  the  slopes 
of  these  five  curves,  except  the  lowest  one  corresponding  to  =  0  dB,  are  less  than  one. 

Recall  the  parabolic  slice  characterizations  for  both  the  equal-arrival  special  case  and  the 
general  case  over  the  vertical  axis  of  the  isosceles-triangle- like  region.  For  the  equal-arrival  special 
case,  the  slope  of  A2  with  respect  to  77  is  00  at  the  vertex  and  is  i  at  the  baseline.  For  the  general 
case  with  >  10  dB,  the  slope  is  about  1  throughout  the  entire  region. 

The  straight  lines  with  slope  of  -f  1  for  ^  >  10  dB  gives  another  indication  of  the  additive  or 
multiplicative  dependency  of  the  total  stress  on  the  two  contributing  factors  and  77  mentioned 
above. 

Mesh  Surface  and  Contour  Plots  of  Signal  Eigenvalue  Ratio  j^(dB)  over  Phase- 
Independent  Variable  T7(dB)  and  Power  Ratio  ^(dB)  for  £=0 

Figures  5.31  and  5.32  for  £  =  0  are  the  mesh  surface  and  contour  plots  of  ^-(dB)  vs  77(dB)=— 40: 1 :0 
and  £(dB)=0:l:30. 

The  previous  two  plots  just  discussed  are  the  slices  resulting  from  vertical  cutting  planes 
intersecting  the  mesh  surface  as  displayed  here.  The  contour  plot  is  that  resulting  from  using 
horizontal  cutting  planes.  The  lack  of  curvature  or  the  straight  contour  lines  for  ^  >  10  dB 
and  the  essential  equal  spacing  of  the  contour  lines  reinforces  the  earlier  discussions.  We  have 
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discussed  in  subsection  3.5.4  that  for  large  £1,  the  contours  of  the  eigenvalue  ratio  =  /  for 
£  =  0  as  functions  of  ^  and  77  are  45°  parallel  lines  with  all  three  variables  ^  =  /,  and  77  in 

dB  scale.  All  these  suggest  that  a  3-d  map  with  power  ratio  ^  as  third  coordinate  is  not  needed. 
The  following  plots  are  the  ”77  in  direct  scale”  versions  of  the  above  ”77  in  dB-scale”  plots. 

5.4.2  With  Phase-Independent  Variable  77  in  Direct  Scale 

Signal  Eigenvalue  Ratio  j^(dB)  vs  Phase-Independent  Variable  77  Curves  at  Seven 
Constant  Power  Ratio  ^(dB)  Values  for  £  =  0 

Figure  5.33  for  £  =  0  has  seven  curves  of  ^(dB)  vs  77=0.05:0.05:1.00  at  ^-(dB)=0,  5,  10,  15,  20, 
25,30. 

Signal  Eigenvalue  Ratio  y^(dB)  vs  Power  Ratio  ^(dB)  at  Five  Constant  Phase- 
Independent  Variable  77  Values  for  £  =  0 

Figure  5.34  for  £  =  0  has  five  curves  of  jh-(dB)  vs  ^(dB)=0:l:30  at  77=0.05,  0.25,  0.50,  0.75,  1.00. 

Mesh  Surface  and  Contour  Plots  of  Signal  Eigenvalue  Ratio  ^(dB)  over  Phase- 
Independent  Variable  77  and  Power  Ratio  |^(dB)  for  £=0 

Figures  5.35  and  5.36  for  £  =  0  are  the  mesh  surface  and  contour  plots  of  ^(dB)  vs  77=0.05:0.05:1.00 
and  fi(dB)=0:l:30. 

We  find  that  the  earlier  ”77  in  dB-scale”  plots  convey  information  more  effectively  and  com¬ 
pactly  than  these  ”77  in  direct-scale”  ones.  But  the  latter  ones  can  be  considered  as  magnified  look 
of  the  earlier  ones  of  the  relatively  benign  region  near  the  apex  with  77  >  -ts  10  dB  or  >  0.10. 
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Chapter  6 

Some  Cadzow’s  and  Ottersten’s 
Scenarios 


In  this  chapter,  we  will  assess  some  scenarios  used  by  Cadzow[8]  and  Ottersten[37]  in  contrasting 
their  results  reported  on  performance  comparison  between  different  ds.  We  first  summarize 

some  approximate  results  derived  in  the  preceding  chapters  that  will  be  used  in  later  example 
scenario  discussions.  These  easy- to- remember  formulas  enable  one  to  exercise  back-of-envelope 
assessments  of  scenarios. 


6.1  A  Concise  Review  of  Approximate  Formulas  Derived  Early 

Among  the  three  parameters  (^,  77,  appearing  in  the  eigenvalues’  expressions,  77  is  the  most 
important  one.  The  many-to-one  mapping  from  (p,  <f>)  to  (^, 77),  or  from  the  Cartesian  product 
of  the  temporal  and  spatial  correlation  unit-disks  onto  the  isosceles-like  region,  is  defined  in 
section  3.1  as 

7  =  (1  -  l<£l2)(i  -  H2)> 

£  =  Re(pd>*)  =  \p\\4>\  cos(argp  —  arg(£). 

When  the  scenario  is  stressful,  7?  is  small  and  the  corresponding  annuli  with  unit  outer  radius  are 
thin.  The  area  of  the  temporal  correlation  p  annulus  can  be  approximated  by  the  following  (aside 
from  a  factor  of  tt): 

l-|p|2  =  (l  +  |p|)(l-|p|)«2(l  -  |p|). 

For  two  closely  spaced  plane  waves  impinging  on  a  uniform  line  array  (ULA),  we  have  the  ap¬ 
proximate  formula 

V  ~  jM/rac.hw)2- 
or 

*?(in  dB)  «  5  +  2  x  9frac_bw( in  dB), 

as  discussed  in  section  4.2.  The  formula  is  at  least  90%  accurate  regardless  of  the  number  of 
sensor  elements. 

For  —  0  dB  and  rj  <  -10  dB,  we  have  the  approximate  formula  for  £  =  0 

^-(in  dB)  «  -r?(in  dB)  +  6, 

A2 
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as  derived  in  section  3.3.  For  these  highly  correlated  equipower  arrivals,  the  effect  on  the  eigenvalue 
ratio  from  changing  the  phase  of  the  source  correlation  from  90°  to  0°  can  be  as  high  as  a  6-dB 
increase,  as  derived  in  section  6.3. 

For  ^  >  10  dB,  we  have  the  approximate  formula  for  £  =  0 

^-(in  dB)  «  (— )(in  dB)  -  7?(in  dB), 

A2  7T2 

as  derived  in  section  3.3.  The  effect  on  the  eigenvalue  ratio  from  changing  the  phase  of  the  source 
correlation  from  90°  to  0°  is  negligible  for  ^  >  15  dB  as  discussed  in  section  5.3. 

6.2  Temporally  Uncorrelated  Arrivals 

The  most  stringent  scenario  used  by  Cadzow  for  two  nearly  broadside  uncorrelated  arrivals  im¬ 
pinging  on  a  circular  arc  array  of  eight  sensors  equally  spaced  at  3°  is  now  discussed.  For  /o=15 
MHz,  a  300-m  diameter  is  15  wavelengths.  As  the  circular  arc  is  nearly  linear,  the  spatial  correla¬ 
tion  of  the  steering  vectors  of  the  arrivals  for  the  exact  circular  arc  array  is  estimated  by  that  of 
an  approximating  ULA.  A  24°-arc1  with  a  7.5- wavelength  radius  presents  an  effective  aperture  of 

7.5  wavelength  x24°  X  =  x  wavelength 

for  broadside  DOAs.  The  4-dB  down  Raleigh  beamwidth  4^dB.down.bwJn.0  is  18.2°.  The  closest 
spacing  between  two  arrivals  used  by  Cadzow  is  3°  which  corresponds  to  a  fractional  beamwidth 
of  ^jrac-bw  =  0.16  for  an  approximating  ULA. 

For  a  ULA  with  eight  sensors  and  fractional  beamwidth  9frac^bw  =  0.16,  we  read  off  from 
figure  4.3  the  value  of  A2  =  1  -  |0|  =  -14  dB  or  77  -  -11  dB  as  we  noted  A2  ~  |  in  section  4.2  for 
uncorrelated  and  equal  strength  arrivals.  The  accuracy  of  reading  off  curves  is  probably  no  better 
than  using  the  approximate  formula  for  closely  spaced  arrivals  which  yields  |(7r  x  0.16)2  « 
or  -11  dB.  Then  we  read  from  the  curve  for  ^  =  0  dB  in  figure  5.33  with  77  in  direct  scale 
the  value  of  ^  =  16  dB.  This  can  also  be  obtained  from  the  approximate  formula  which  is 
valid  for  ^  =  0  dB  and  77  less  than  -10  dB.  So,  the  optimum  weight  ratio  is  like  40:1.  The 
tradeoffs  conducted  by  Cadzow  were  mainly  between  his  (unweighted)  algorithm  and  Det-MLE. 
The  unweighted  algorithm  was  outperformed  by  Det-ML,  but  not  by  much. 

The  ^  ratios  used  by  Ottersten  for  his  ULA  scenarios  are  typically  in  the  area  of  30  dB,  i.e., 
1000:1,  or  higher.  For  the  two  equal-SNR  uncorrelated  arrivals  within  a  Raleigh  beamwidth,  the 
most  stringent  beamwidth  separation  is  6]TacJbw  —  1/40  =  0.025.  This  is  for  the  three-element 
half-wavelength  spacing  case  with  only  1°  separation.  For  a  ULA  with  three  sensors  and  fractional 
beamwidth  0fracJ> w  —  0.025,  we  read  off  from  figure  4.1  the  value  of  A2  =  1  —  |(^|  =  —30  dB,  i.e., 
77  =  —27  dB.  This  can  also  be  obtained  from  the  approximate  formula  for  closely  spaced  arrivals 
as  f(7r/40)2  «  ^ .  Then  we  read  from  the  curve  for  =  0  dB  in  figure  5.29  with  77  in  dB  scale 
the  value  of  ^  =  32  dB,  or  1600:1,  which  is  in  agreement  with  the  approximate  formula  for  ^  =  0 
dB  and  77  less  than  0.10  (or  -10  dB).  Ottersten’s  results  indicate  that  the  degradation  from  using 
unit  weights  is  a  five-time  increase  in  the  DOA  standard  deviation  for  the  four  sensors  case  and 
a  three-time  increase  for  10  sensors  case. 

The  next  uncorrelated  scenario  used  by  Ottersten  is  for  two  unequal-SNR  arrivals  at  30  and 
0  dB  and  spaced  much  farther  apart  at  0°  and  15°  from  the  broadside.  For  the  most  stringent 

‘Recall  the  discussions  from  section  4.1,  the  effective  aperture  of  the  approximating  ULA  should  be  a  24°  arc 
length  instead  of  a  21°  one. 
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case  of  a  three-element  half-wavelength  spacing  ULA,  this  is  about  40%  fractional  beamwidth. 
From  figure  4.3,  we  find  from  77  =  1  —  \<f>\2  vs  fractional  beamwidth  for  the  three-element  case  the 
value  of  Tj  —  — 3  dB,  roughly.  This  can  also  be  obtained  from  the  approximate  formula  for  closely 
spaced  arrivals  as  |(7r  X  0.4)2  «  0.48.  Then,  we  read  from  the  curve  for  ^  =  0  dB  in  figure  5.33 
with  77  in  direct  scale  the  value  of  ^  =  33  dB.  This  can  also  be  obtained  from  the  approximate 
formula  which  is  valid  for  ^  >  10  dB.  So,  the  optimum  weight  ratio  is  about  2000:1. 

In  addition  to  the  higher  DOA  standard  deviation  for  using  the  unit  weights,  the  much  more 
offending  and  eye-catching  phenomenon  [38,  page  68]  is  the  roller-coaster  behavior  of  the  standard 
deviation  curve  as  a  function  of  the  number  of  sensors.  That  is,  the  performance  can  degrade 
by  adding  more  sensors.  This  has  been  found  by  Ottersten  [38,  page  69]  as  highly  correlated 
with  the  sidelobe  behavior  of  the  ULA.  Unlike  the  previous  example  where  the  two  arrivals  stay 
within  the  mainlobe,  here  each  arrival  gets  farther  and  farther  into  the  sidelobes  of  the  other 
as  more  sensors  are  used.  Our  plot  of  the  normalized  inner  product  of  the  two  steering  vectors 
agree  with  Ottersten ’s  plot,  but  the  ^  ratio  curve  has  the  maxima  of  34.5  dB  and  exhibits  only 
0.5  dB  ripples.  Evidently,  the  asymptotic  standard  deviation  expression  is  closely,  if  not  exactly, 
proportional  to  this  normalized  inner  product.  This  is  beyond  the  scope  of  our  approach,  where 
only  noise-free  simplification  is  used. 

As  this  defect  of  the  unit  weight  has  the  close  relation  to  the  sidelobe  phenomena,  it  probably 
will  not  be  prominent  for  a  random  array.  However,  the  drawback  of  not  using  appropriate  weight¬ 
ing  will  be  always  present  and  manifested  by  the  deterioration  of  the  DOA  standard  deviation, 
especially  when  high  ^  contrast  is  there. 

6.3  Temporally  Highly  Correlated  Arrivals 

Ottersten  provided  two  additional  scenarios  using  equal- SNR  partially  but  highly  correlated  ar¬ 
rivals.  This  is  contrasted  to  the  either  0%-  or  100%-correlated  and  equal-SNR  arrivals  typically 
used  by  Cadzow.  Here  Ottersten’s  results  are  displayed  as  standard  deviation  versus  SNR.  We 
need  to  know  more  about  the  multipath  correlation  beyond  the  magnitude,  especially  for  the  case 
of  equal  SNRs  where  the  variable  £  involving  the  phase  term  has  its  most  effect.  During  the 
discussion  of  these  scenarios,  we  will  first  tentatively  assume  that  the  phase  difference  between 
the  two  equal-SNR  multipaths  as  90°  and  at  the  end  adjust  the  eigenvalue  ratio  as  we  move  to 
the  0°  correlation.  For  these  highly  correlated  equipower  arrivals,  we  have  shown  in  section  5.3 
that  the  effect  on  the  eigenvalue  ratio  from  changing  the  phase  of  the  source  correlation  from  90° 
to  0°  can  be  as  high  as  a  four-time  or  6-dB  increase. 

Partial  correlation  of  99%  yields  a  partial  contribution  to  77  of  -17  dB,  as  1  -  j0.99|2  % 
2(1  -  |0.99|)  =  0.02.  The  separation  of  two  DOAs  of  7°  for  a  half- wavelength  spacing  ULA 
corresponds  roughly  to  0jrac_bw  =  0.2  beamwidth,  which  contributes  about  -9  dB  stress  as  read 
off  from  figure  4.3.  This  can  also  be  obtained  from  the  approximate  formula  for  closely  spaced 
arrivals  as  |(7r  x  0.2)2  «  0.12.  Therefore,  we  have  77  =  -26  dB  because  77  is  exactly  multiplicative, 
or  additive  in  dB,  in  its  two  contributing  factors,  i.e.,  the  multipath  correlation  part  and  the 
closeness  of  the  two  steering  vectors  of  the  DOAs.  So,  the  eigenvalue  contrast  here  for  the  noise- 
free  scenario  can  be  read  from  the  curve  for  =  0  dB  in  figure  5.29  using  77  =  -26  dB  as  ^  =  32 
dB  or  1600:1.  This  is  in  agreement  with  the  approximate  formula  valid  for  ^  =  0  dB  and  77 
less  than  -10  dB.  After  the  full  4-dB  adjustment  needed  for  the  0°  temporal  correlation,  we  have 
^  =  36  dB  or  4000:1.  For  moderate  SNRs  such  as  those  about  0  dB  ones,  the  deterioration  of 
the  Det-ML  with  respect  to  WSF  starts  to  be  noticeable  as  compared  to  the  lack  of  distinction 
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between  the  two  at  higher  SNRs.  This  will  be  elaborated  in  the  next  paragraph.  The  degradation 
of  the  unit  weighting  from  that  of  the  optimal  one  is  more  than  one  order  of  magnitude. 

Partial  correlation  of  90%  yields  a  partial  contribution  to  r?  of -7  dB,  as  1  -  |0.90|2  ~  2(1  - 
|0.90|)  =  0.2.  The  40°-wide  separation  of  two  DOAs  corresponds  roughly  to  1.3  beamwidth  for 
a  four-sensor  ULA  which  contributes  about  0-dB  stress  as  shown  from  figure  4.7.  So,  the  stress 
is  almost  entirely  from  the  partial  correlation  in  the  noise-free  simplification.  Therefore,  we  have 
T)  —  -7  dB  or  0.2.  This  is  in  agreement  with  the  approximate  formula  valid  foT  =  0  dB  and 
t ;  less  than  -10  dB.  The  eigenvalue  contrast  here  for  the  noise-free  scenario  can  be  read  from  the 
curve  for  ^  =  0  dB  in  figure  5.33  with  tj  in  direct  scale  using  t)  =  —7  dB  as  ^  =  13  dB  or  200:1. 
Ottersten  used  fairly  low  SNR  but  a  large  number  of  snapshots  to  distinguish  the  optimal  weights 
were  used  from  the  more  intuitive  MLE  weights.  For  large  SNRs,  these  two  different  choices  of 
weights  do  not  have  significantly  different  values.  Therefore,  this  choice  may  not  be  important  in 
the  HFDF  scenarios.  But  as  SNRs  decrease,  the  magnitude  of  these  two  different  sets  of  weights 
depart,  and  the  superiority  of  the  optimal  WSF  algorithm  over  the  Det-ML  gets  more  pronounced. 
Any  preference  of  the  optimal  weights  over  the  more  intuitive  MLE  weights  in  the  passive  sonar 
context  is  also  questionable:  to  estimate  the  location  of  threshold  signals,  the  number  of  signals 
must  be  determined.  None  of  these  highly  parametric  direction-finding  algorithms  is  suitable  for 
threshold  signals. 

6.4  Temporally  100%-Correlated  Arrivals 

So  far,  we  have  contrasted  the  Monte  Carlo  simulation  comparison  results  of  some  of  Cadzow’s 
0%-correlated  scenarios  and  the  comparison  results  of  both  the  Monte  Carlo  simulation  and 
the  asymptotic  performance  prediction  of  some  Ottersten ’s  uncorrelated  and  partially  correlated 
scenarios.  We  found  it  useful  to  assess  the  severity  of  their  scenarios  in  terms  of  the  noise-free 
eigenvalue  distribution  for  the  two  discrete  arrivals. 

However,  such  assessment  fails  when  the  small  eigenvalue  for  the  noise-free  simplification 
is  zero,  such  as  the  case  for  the  100%-correlated  arrival  case  which  was  discussed  by  Cadzow. 
Secondly,  there  is  some  qualitative  discrepancy  between  the  asymptotic  and  simulation  results  for 
the  tradeoffs  for  high-SNR-scenario  cases  for  the  100%-correlation  cases  between  the  Det-ML  and 
WSF  methods.  The  asymptotic  results  cannot  differentiate  the  performances  of  these  two  methods 
by  just  judging  from  the  weight  expressions.  Ottersten  did  not  provide  any  asymptotic  and 
simulation  results  for  the  100%-correlation  cases.  Whereas  Cadzow’s  simulation  results  for  some 
100%-correlated  arrival  cases  indicated  that  even  the  unit  weighting  signal  eigenvector  method 
outperform  the  Det-ML  method.  For  the  single  dominant  eigenvector  case,  there  is  no  weighting 
involved  at  all,  consequently  there  should  be  no  difference  between  Cadzow’s  signal  eigenvector 
method  and  the  WSF  method. 

As  the  asymptotic  Det-ML  for  the  high-SNR  case  approaches  the  CRLB,  it  is  puzzling  to  see 
how  this  bound  for  an  unbiased  estimator  can  be  exceeded.  Ottersten  believes  that  apparently 
this  is  because  of  the  longer  time  it  takes  for  the  finite  sample  size  bias  to  diminish  in  Det-ML  than 
WSF.  Asked  to  provide  some  plausible  explanations  for  the  reason  that  the  signal  eigenvector  type 
methods  outperform  the  Det-ML  one,  Ottersten  thinks  that  this  is  because  of  the  rank-reduction 
advantage  provided  in  using  only  the  dominant  eigenvectors.  While  the  original  rank-reduction 
idea  was  suggested  by  Scharf  and  Tufts  [47]  and  Scharf  [48]  in  some  different  context,  they  did 
point  out  the  potential  applications  to  DF. 

We  can  see  readily  that  Cadzow’s  results  fitting  two  100%-correlated  signals  with  a  single 
eigenvector  was  not  bad  at  all.  f  urthermore,  when  only  one  eigenvector  is  involved,  the  weighting 
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of  eigenvectors  is  irrelevant,  i.e.,  the  unweighted  signal  eigenvector  and  WSF  methods  are  identical. 
On  the  other  hand,  the  results  from  using  nonlinear  least  squares  fit  to  the  two  unweighted 
eigenvectors  as  reported  in  the  two  scenarios  of  the  90%-  or  99%-correlated  DOAs  by  Ottersten  was 
disastrous.  Ottersten  pointed  out  that  in  these  cases,  any  consistent  detection  scheme  will  detect 
two  arrivals  given  enough  samples.  He  also  pointed  that  there  is  a  discontinuity  in  performance 
at  100%  correlation  for  the  unweighted  signal  eigenvector  method,  whereas  WSF  does  not. 

When  a  ’’dominant  eigenvalue”  becomes  difficult  to  be  distinguished  from  noise  level  ones, 
one  would  intuitively  choose  not  to  use  it.  This  is  based  on  the  conjecture  that  retaining  this 
marginally  signal  eigenvector  and  using  the  nonlinear  least  squares  fit  to  it  with  an  additional 
linear  combination  of  steering  vectors  corresponding  to  the  unknown  DOAs  probably  will  do  more 
harm  than  help.  Scharf  and  Tufts  gave  such  practices  more  legitimate  reasoning.  Consequently, 
the  problem  of  deciding  the  optimal  number  of  signal  eigenvectors  to  be  used  in  practice  may  not 
be  a  simple  one. 


Chapter  7 

Conclusions 


After  describing  the  scope  of  the  report,  chapter  1  starts  with  a  survey  of  NB  DF  algorithms 
for  general  sensor  arrays  in  element  space.  This  includes  a  taxonomy  table,  characterization  and 
asymptotic  performance  of  the  main  algorithm  categories,  and  the  initialization  procedure  which 
is  the  critical  item  of  most  high-resolution  algorithms.  Then,  formulating  guidelines  desired  for  DF 
scenario  design  are  developed  via  the  following  items:  too  many  parameters  for  a  two-source  DF 
scenario  for  a  given  sensor  array,  simple  map(s)  to  show  relative  positions  of  scenarios,  yardsticks 
for  measuring  stress  of  scenarios  to  direction-finding  algorithms,  and  dominant  eigenvalues  chosen 
as  the  key  items.  Next,  we  discuss  eigenvalues  for  the  noise-free  covariance  matrix  of  two  discrete 
sources  through  the  following  issues:  treatment  of  this  subject  in  Hudson’s  text,  the  original  6- 
dimensional  parameter  space  for  a  given  sensor  array  under  the  noise-free  assumption,  compact 
map(s)  -  another  many-to-one  reduction,  and  3-d  map  with  power  ratio  ^  as  third  coordinate. 

Chapter  2  addresses  the  notion  of  stress  to  direction-finding  algorithms  and  scenarios.  First, 
the  three  main  sources  for  the  signal  eigenvalue  spread  or  causes  for  small  signal  eigenvalues 
are  reviewed:  high  spatial  and/or  temporal  correlation  and/or  high  contrast  of  signal  strength 
between  the  two  arrivals.  Next,  the  eigenvalues  in  weighting  expressions  are  discussed  in  three 
categories:  in  noise  subspace  methods,  in  parametric  signal  subspace  fitting  methods,  and  in  non- 
parametric  subspace  adaptive  beamforming.  Then,  the  drawbacks  of  the  oversimplifying  stress- 
measuring  yardstick  using  condition  numbers  are  mentioned:  infinite  condition  number  for  all 
100%-correlated  two  arrivals,  equal  condition  numbers  corresponding  to  different  CRLBs,  CRLB 
not  appealing  to  intuitive  expectation  sometimes,  and  the  redeeming  value  for  condition  number. 
At  the  end,  we  explain  the  reasons  why  special  treatment  is  given  for  equal-power  arrivals  and 
the  complementary  value  of  direct-scale  and  dB-scale  plots. 

The  analytical  chapter  3  begins  by  reviewing  the  expressions  of  the  eigenvalues  of  the  quadratic 
characteristic  equations  of  the  non-Hermitian  product  of  the  temporal-  and  spatial- correlation 
matrix.  The  eigenvalues  are  normalized  with  respect  to  the  product  of  the  number  of  sensors 
and  the  sensor  level  power  of  the  weaker  source  ir?.  The  normalized  large  and  small  eigenvalues 
expressions  A1  and  A2  are  reduced  to  forms  depending  only  on  the  real  triplet  (£,  r?,  JJ-).  Here,  jp- 
is  tb*>  power  ratio  between  the  strong  and  weak  sources  at  the  sensor  level.  The  real  pair  (£,  t?) 
are  defined  in  terms  of  the  normalized  temporal-  and  spatial-coefficients  p  and  </>  respectively  with 
each  constrained  to  a  unit-disk  through 

^(l-|^|2)(l-|p|2), 
i  =  Re(p</>*)  =  |p||d>j  cos(argp  -  argd>). 
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While  the  discussions  of  eigenvalues  in  this  report  are  applicable  to  generic  steering  vectors,  we 
call  77  the  phase-independent  variable  and  £  the  phase-dependent  variable  for  convenience  even 
though  such  names  were  motivated  by  the  special  case  of  plane  waves  impinging  onto  pairwise 
symmetric  arrays.  The  phase-independent  variable  77  depends  only  on  the  magnitude  of  p  and  4>. 
The  phase-dependent  variable  £  depends  aioO  on  the  angular  positions  of  p  and  4>.  The  mapping 
from  the  complex  pair  (p,  <f>)  to  the  real  pair  (£,77)  is  many-to-one.  The  range  of  this  mapping  is 
an  isosceles-right-triangle-like  region  bounded  on  its  left  and  right  by  two  symmetric  parabolas 
meeting  each  other  at  90°  and  down  below  by  a  straight  baseline.  We  characterize  this  region 
thoroughly  with  respect  to  the  temporal  and  spatial  coefficients  of  direction-finding  scenarios. 
Most  of  this  report  is  devoted  to  characterizing  the  behavior  of  the  eigenvalues  and  their  ratio 
as  functions  of  the  real  triplet  (£,77,^).  For  convenience,  we  jointly  summarize  chapter  3, 
which  is  mainly  on  the  eigenvalue  analysis,  and  chapter  5,  which  is  mainly  on  the  eigenvalue 
ratio’s  display.  After  the  expressions  for  the  eigenvalue  and  their  ratio  ^  for  the  special  cases 
at  the  apex,  the  baseline,  and  the  vertical  axis  were  discussed,  we  thoroughly  characterize  the 
special  case  of  equal-strength  arrivals  and  then  the  general  case. 

The  special  case  of  equal  strength  =  0-dB  arrivals  share  the  following  common  features  of 
general  -j-  power  ratios:  straight  line  contours  for  constant  eigenvalue  A,  parabolic  A  slices  for 
constant  phase-dependent  variable  £,  hyperbolic  A  slices  for  constant  phase-independent  variable 
77,  and  parabolic  contours  for  eigenvalue  ratio 

The  special  case  of  equal  strength  —  0-dB  arrivals  merits  unique  attention  because  of  the 
following.  There  is  an  extra  6  dB  for  signal  eigenvalue  ratio  ^  over  the  vertical  axis  £  =  0, 
accounting  for  the  presence  of  an  infinite  slope  of  A  with  respect  to  77  at  the  apex  for  the  equal- 
power  arrival  case.  This  is  reflected  in  the  following  two  approximations.  When  ^  >  10  dB,  we 
have 

~{m  dB)  «  (  — )(in  dB)  -  r?(in  dB). 

A2  7T2 

When  ^  «  0  dB  and  77  is  less  than  0.10  (or  -10  dB),  we  have 

Y^-(in  dB)  ss  -7?(in  dB)  +  6. 

A2 

Furthermore,  for  the  special  case  of  equal  strength  ^  =  0  dB  arrivals,  the  contours  of  the 
eigenvalue  ratio  ^  =  /  form  a  one-parameter  family  of  parabolas  having  their  vertices  colocated 
at  (£,77)  =  (—1,0),  i.e.,  at  the  lower  left  corner  of  the  isosceles-like  triangle  and  their  common 
tangent  the  baseline  77  =  0.  This  family  include  the  left  boundary  of  the  isosceles-like  triangle 
and  the  baseline.  Because  the  contours  of  the  eigenvalue  ratio  ^  =  /  pass  through  the  lower  left 
corner  of  the  isosceles-like  triangle,  this  corner  is  a  point  of  discontinuity  for  the  eigenvalue  ratio 
for  the  special  case  of  equal  strength  ^  =  0  dB.  For  example,  the  two  eigenvalues  are  equal 
over  the  left  parabolic  boundary  of  the  isosceles-right-triangle-like  region.  Furthermore,  the  two 
equally  dominant  eigenvalues  can  go  down  to  zero  jointly  so  that  the  condition  number  stays  as 
unity  throughout.  That  is,  the  condition  numbers  are  the  lowest  possible  there.  However,  both 
eigenvalues  tend  to  zero  as  77  tends  to  zero,  i.e.,  the  baseline  of  the  region. 

As  power  ratio  ^  increases  from  unit  value,  the  variation  of  the  large  eigenvalue  A]  is  es¬ 
sentially  along  the  direction  of  the  horizontal  phase-dependent  variable  £,  and  that  the  small 
eigenvalue  A2  essentially  along  the  direction  of  the  vertical  phase-independent  variable  77. 

There  is  a  four-time  or  6-dB  increase  on  the  eigenvalue  ratio  ^  for  highly  correlated  equipower 
arrivals,  both  temporally  and  spatially,  from  changing  the  phase  of  the  source  correlation  or  the 
angle  difference  between  the  two  unit-disk  vectors  from  90°  to  0°. 
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The  smaller  eigenvalue  A2  is  independent  of  both  the  horizontal  phase-dependent,  variable  £ 
and  ^  near  the  baseline  for  large  but  is  a  strong  function  of  the  vertical  phase-independent 
variable  77.  The  effect  of  (p-  is  essentially  only  felt  by  the  large  eigenvalue  which  is  not  sensitive 
to  whether  the  scenario’s  77)  coordinate  is  close  to  the  baseline,  i.e.,  whether  the  arrivals  are 
correlated  or  close  to  each  other  as  far  as  its  order  of  magnitude  is  concerned. 

We  also  characterize  Speiser’s  eigenvalue  bounds  specialized  to  the  two-source  scenario.  These 
bound  expressions  involve  only  ihe  absolute  magnitude  of  the  two  phasors  on  the  two  unit-disks, 
i.e.,  only  the  radius  but  not  the  argument  of  the  phasors.  Interestingly  enough,  his  two  bounds 
for  the  equipower  case  can  be  identified  as  the  eigenvalue  values  at  the  two  extreme  ends  of  the 
permissible  phase-dependent  variable  ^-interval  on  the  constant  phase-independent  77  horizontal 
line  of  the  isosceles  region  with  all  expressions  calculated  from  the  magnitudes  of  p  and  <f>  only.  For 
given  77,  the  permissible  ^-interval  attains  its  maximum  width  when  the  ^-interval  extends  from 
the  left  to  the  right  boundary  parabola.  The  length  of  such  maximum  ^-interval  increases  with 
decreasing  77,  i.e.,  when  the  baseline  is  approached.  The  bounds  for  the  small  eigenvalue  may  be 
loose  when  there  is  considerable  variation  in  the  order  of  magnitude  of  the  small  eigenvalue  in  the 
phase-dependent  variable  f.  For  large  tp-,  the  smaller  eigenvalue  A2  becomes  more  independent 
of  the  horizontal  phase-dependent  variable  £,  the  bounds  will  get  tighter.  That  is,  there  is  less 
information  loss  by  neglecting  the  phase-dependent  variable  £. 

Even  though  we  assess  only  those  algorithms  applicable  to  gener?1  sensor  array  configurations, 
we  still  need  to  understand  the  results  about  plane  waves  impinging  on  a  ULA,  especially  those 
simulation  results  which  could  be  expensive  to  duplicate  for  general  arrays.  Chapter  4  discusses 
fractional  beamwidth  measurement  for  arrival  separation  and  the  smaller  eigenvalue  A2’s  negligible 
dependence  on  the  number  of  elements  N  after  using  0fr acj,w  for  the  equipower  uncorrelated  case. 

In  chapter  6,  we  apply  all  the  previously  developed  results  to  assess  some  scenarios  used  by 
Cadzow  and  Ottersten.  For  temporally  uncorrelated  arrivals,  we  found  that  the  ^  ratios  used  by 
Ottersten  for  his  ULA  scenarios  are  typically  in  the  area  of  30  dB,  i.e.,  1000:1,  or  more,  whereas 
the  mill-cross  example  used  by  Cadzow  need  an  optimum  weight  ratio  like  40:1.  This  explains 
why  different  conclusions  were  drawn  by  them.  The  tradeoffs  conducted  by  Cadzow  were  mainly 
between  his  (unweighted)  algorithm  and  Det-ML.  The  unweighted  algoiillnii  was  outperformed 
by  Det-ML  but  not  by  much.  Ottersten’s  results  indicate  that  the  degradation  from  using  unit 
weights  is  a  5-time  increase  in  the  DOA  standard  deviation  for  a  4-sensor  case  and  a  3-time 
increase  for  a  10-sensor  case.  However,  in  addition  to  the  higher  DOA  standard  deviation  for 
using  the  unity  weights,  the  much  more  offending  and  eye-catching  phenomenon  is  the  roller¬ 
coaster  behavior  of  the  standard  deviation  curve  as  a  function  of  the  number  of  sensors.  That 
is,  the  performance  can  degrade  by  throwing  more  sensors.  This  has  been  found  by  Ottersten  as 
highly  correlated  with  the  sidelobe  behavior  of  the  ULA. 

We  also  commented  on  Ottersten’s  scenarios  of  temporally  highly  correlated  arrivals  and 
C’adzow’s  scenarios  of  temporally  100%-correlated  arrivals.  For  the  latter  case,  the  eigenvalue  ratio 
^  =  00  and  the  attempt  to  assess  the  severity  of  their  scenarios  via  the  noise-free  eigenvalue  ratio 
for  the  two  discrete  arrivals  fails.  There  is  some  qualitative  discrepancy  between  the  asymptotic 
and  simulation  results  for  the  tradeoffs  for  high-SNR  scenario  cases  for  the  100%  correlation 
rases  between  the  Det-ML  and  WSF  methods.  The  asymptotic  results  cannot  differentiate  the 
performances  of  these  two  methods  just  judging  from  the  weight  expressions.  Ottersten  did  not 
provide  any  asymptotic  and  simulation  results  for  the  100%  correlation  cases.  Whereas  Cadzow’s 
simulation  results  for  sonm  100%-correlated  arrival  cases  indicated  that  even  the  unity  weighting 
signal  eigenvector  method  outperform  the  Det-ML.  Ottersten  conjectured  that  the  rank-reduction 
idea  by  Scharf  and  Tufts  may  answer  such  different 
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Appendix  A 

Hyperbola  Derivation  and 
Characterization 


In  the  general  case  discussion,  we  asserted  without  proof  that  for  constant  77  slices  the  characteristic 
equation  can  be  rearranged  to  that  of  a  hyperbola  in  A  and  £  as 

A2  -  ,/^2£A  -  (1  +  — )A  +  —t)  =  0. 

V  n2  7T2  *2 

In  this  appendix,  we  will  show  that  for  each  ^  value,  this  equation  describes  a  one-parameter  77 
family  of  hyperbolas  in  A  and  £,  sharing  common  asymptotes. 

Recall  that  in  the  discussion  of  the  special  case  of  77  =  0,  i.e.,  at  the  baseline,  we  have  the 
following  relations: 

A2  =  0, 

A,  =  —  +  1  +  ^  2£. 


When  £  changes  from  -1  to  0  to  + 1,  Ai  varies  linearly  in  (  with  the  slope  equal  to  from 


(^_1)!to(ji  +  1)lo(^  +  i). 

These  two  straight  lines  for  Aj  and  A2  will  turn  out  to  be  the  common  asymptotes  of  each  family 
labeled  with  the  particular  ^  value. 

Following  the  usual  treatment  for  reducing  a  general  second-order  plane  curve  to  its  canonical 
form,  we  make  the  following  coordinate  transform 


cos  9  -  sin  9 
sin#  cos  9 


E  cos  9  -  A  sin  9  +  h 
E  sin  9  +  A  cos  9  +  k 


Substituting  into  the  original  equation,  we  have 


(^.sin  9  +  A  cos  9  +  kj2  -  «/— 2(Ecos0  -  A  sin  9  +  /i)(E  sin  9  +  A  cos  9  +  k) 

V  *2 

-(1  +  —  )(Esin0  +  Acos0  +  k)  +  —  77  =  0. 

7T2  *2 


A-l 


This  can  be  rewritten  into  the  following  form  as 

C \ 2A2  +  C=a— A  +  (7=2— 2  +  C\A  +  Ch—  4-  c  =  0, 


with 


CA 2  =  cos2  0  +  2,/ — sin  0  cos  0, 

CS  A  = 


*  1, 


2sin  0cos0  -  2. /—(cos2  0  -  sin2  0), 

V  *2 
*7 


(7=2  =  sin2  5  -  2,/—  sin  (9  cos 

V  *2 

Ca  =  2k  cos  6  —  2./— (— fcsint?  f  hcosO)  -  (1  H — -)cos0, 

V  ?r2  7r2 

(7=  =  2fcsin  0  —  2« /-T(fc  cos0  + /isin  0)  —  (1  4 — -)sin0, 

V  7r2  ^2 

c  =  A:2  -  2,  f^-hk  -  (1  4-  —  )A:  +  —  7;. 

V  ^"2  ^"2 

The  rotation  of  coordinates  is  obtained  by  setting  C= a  =  0,  i.e., 

tan  26  =  2  ,/^X 


Using  this  result,  we  rearrange  Ca^  and  <7=2  to  see  that  what  we  are  dealing  with  are  indeed 
hyperbolas: 

sin2  28  1 

--  =  1/2  +  77-  77^  >  1, 


<7A2  =  (1  +  cos20)/2  4- 
(7=2  =  (1  -  cos  28)/ 2 


2  cos  20 
sin2  20 


=  1/2  - 


2  cos  20 
1 


<  0. 


2  cos  20  2  cos  20 

The  translation  of  coordinates  is  obtained  by  setting  C\  —  0  and  C-=  =  0,  or 


h< 2  -  2 


-ros'- '•+/?>  =  1  + 7 

7r2  tan  0  V  7t2  7 r2 


7Ti 


fc(2  4-  2y  -T  tan  0)  —  h(2^J =  1  +  ~- 

Therefore  k  =  0,  as  tan 0  ^  We  also  get  h  =  -i(l  +  This  means  that  the  new 

coordinate  origin  is  at  (-^(1  4-  i.e.,  the  intersection  of  the  two  straight  lines  for  A! 

and  A2  introduced  during  the  discussion  of  the  special  case  of  tj  =  0,  i.e  ,  at  the  baseline. 

The  remaining  constant,  c,  is  ^rj  or  ( \  tan  20)2t]. 

The  hyperbola  in  the  new  coordinate  system  has  the  following  form 

.  1  .  ,._2  -  1 


1  2  cos  20 

The  asymptotes  have  the  slopes  of 


1/2)-  -  ( 


,  „  +  1/2)A2  =  —  t). 

2  cos  20  7r2  ' 


-C= 2 

Q2 


'  1  -  cos  20 
1  +  cos  20 


=  ±  tan  0, 


A-2 


as  expected. 

The  vertex  distance  from  the  new  coordinate  origin  is  Ev  obtained  from  setting  A  =  0  as 


_  2  A  i  2  2  cos  2  9  sin  20  1  1 

(-r)  =  (-  tan20)2  ?/- - —  =  - -z  (-tan  25)?/  =  - — jzV 

2  l  —  cos  20  1  —  cos  20  2  1  —  tarr  0 


in  terms  of  0 ,  or 


cos  29  7 r,  (1+  lf^-)  2  Jr  i 

- - — 2  r]  =  - - - : - 2// 

1  -  cos  '28  7r2  i  _  ( i  .}  £v 2  tt2 


in  terms  of 

For  given  0  or  Ev  is  proportional  to  ^/rj.  That  is,  the  rate  of  increase  of  Ev  with  respect 
to  i]  is  oc  at  i)  =  0  and  Ev  at  r;  =  ^  is  already  half  of  Ev  at  rj  =  1.  We  note  that  the  useful 
vertices  of  these  connected  hyperbolas  are  over  the  boundary  of  the  domain  of  interest  only  when 
^  =  1.  Withor*  using  many  data  points  in  the  plotting  and  explicitly  drawing  the  E  a^xis,  these 
aspects  may  not  be  easily  observed  from  the  plots. 
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depending  only  on  a  real  triplet:  phase-dependent  variable  £,  phase-independent  variable  //,  and 
power  ratio  fF  The  pairs  (£.  7/)  are  confined  to  an  isosceles-like  region. 

The  main  contribution  of  this  work  is  a  manageable  presentation  of  a  compact  map  showing 
A 1 .  A 2 .  and  as  functions  of  ?/.  and  ^  over  all  possible  scenarios.  This  enables  one  to  see  the 
relative  positions  among  different  scenarios.  Presented  also  are  some  easy-to-remember  formulas 
that  enable  one  to  exercise  "back-of  envelope"  assessment  of  scenarios. 

The  small  eigenvalue  is  shown  to  diminish  qualitatively  and  quantitatively  for  two-arrival 
scenarios  increasingly  stressed  with  high  temporal  and/or  spatial  correlations.  The  special  case  of 
equal-strength  signal  arrivals  =  1).  also  important  in  low-angle  radar  tracking,  shares  many 
rich  structures  of  general  The  equal-strength  case  also  has  several  additional  unique  features 
for  signal  eigenvalue  ratio  At.  which  is  important  in  direction-finding. 
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